Invasion Percolation Algorithm Research Articles

A computationally efficient queue-based algorithm for simulating volume-controlled drainage under the influence of gravity on volumetric images of porous materials

Simulating non-wetting fluid invasion in volumetric images of porous materials is of broad interest in applications as diverse as electrochemical devices and CO2 sequestration. Among available methods, image-based algorithms offer much lower computational cost compared to direct numerical simulations. Recent work has extended image-based methods to incorporate more physics such as gravity and volume-controlled invasion. The present work combines these two developments to develop an image-based invasion percolation algorithm that incorporates the effect of gravity. Additionally, the presented algorithm was developed using a priority queue algorithm to drastically reduce the computational cost of the simulation. The priority queue-based method was validated against previous image-based methods both with and without the effect of gravity, showing identical results. It was also shown that the new method provides a speedup of 20X over the previous image-based methods. Finally, comparison with experimental results at three Bond numbers showed that the model can predict the real invasion process with a high accuracy with and without gravitational effects.

Pore-scale binary diffusion behavior of Hydrogen-Cushion gas in saline aquifers for underground hydrogen Storage: Optimization of cushion gas type

Cushion gas is required to prevent hydrogen (H2) trapping and interactions between H2 and brine for underground hydrogen storage (UHS) in saline aquifers. However, selecting an appropriate type of cushion gas is not a straightforward task, as it affects the H2 diffusion characteristic in subsurface porous media, thereby impacting purity, recoverability, etc. In this paper, we developed a modified pore network model (PNM) for binary diffusion that incorporates multi-scale diffusion mechanisms to forecast the H2 effective diffusion coefficient (DH2e) and evaluate the impacts of the cushion gas on DH2e. Additionally, the model accounts for the water-lock effects of residual brine, focusing exclusively on the connected pore-throats occupied by the gas phase, which are identified using the invasion percolation algorithm with trapping. Our key findings are as follows: 1) CO2 and N2 demonstrate superior containment effects on H2 compared to CH4 under the specific conditions of the target aquifer, reducing the H2 contamination caused by binary diffusion; 2) The size and topology of porous media exert minimal impacts on the selection of cushion gas, yet they markedly affect the DH2e; 3) The pressure and temperature conditions of the saline aquifer are critical factors in cushion gas selection, with temperature being the more influential parameter. This study provides valuable insights for the implementation of industry-scale UHS in saline aquifers.

Dynamic In Vivo Computation for Learning-Based Nanobiosensing in Time-Varying Biological Landscapes

We have recently proposed a framework of in vivo computation which transforms the early tumor sensing problem into a computational problem. In the framework, a tumor-triggered biological gradient field (BGF) guides swarm-intelligence-assisted targeting process, where externally manipulable and trackable magnetic nanorobots act as computational agents for the optimization procedure. As BGF can be viewed as an objective function which is utilized to define the fitness landscape for the agents, the inherent attributes of BGF are critical to the in vivo computational process. All our previous investigations are based on the hypothesis that the BGF landscape remains time-invariant during the tumor targeting process, which results in a static function optimization problem. However, the properties of internal environment, such as the flow state of body fluid, will naturally lead to time-dependent variation of BGF, which means that the targeting process should be modeled as a dynamic function optimization problem. Based on this consideration, we focus on dynamic in vivo computation by considering different variation patterns of BGF in this paper. Two computational intelligence strategies named “swarm-based learning” and “individual-based learning” are proposed for dealing with the turbulence of the fitness estimation caused by the BGF variation. The in silico experiments and statistical results demonstrate the effectiveness of the proposed strategies. In addition, the above process is conducted in a three-dimensional search space, where the tumor vascular network is generated by an invasion percolation algorithm, which is more realistic compared to the two-dimensional search space in our previous works.

Constraining hydrocarbon migration and accumulation by two-dimensional numerical simulation: Ordovician carbonate reservoirs of the Daniudi Area, Ordos Basin

The various types of reservoir space in carbonate reservoirs make the process of hydrocarbon migration and accumulation very complex. How to quantitatively characterize the characteristics of hydrocarbon migration and accumulation in carbonate reservoirs is one of the key problems for carbonate hydrocarbons exploration and exploitation. Taking the Ordovician in the Daniudi area of Ordos Basin as an example, three types of source rock-reservoir-cap assemblages, namely upper-generation assemblage, lateral variable assemblage, and self-generation-reservoir-cap assemblage, were determined by characterizing the forming elements of carbonate reservoirs. Petromod software is used to build a two-dimensional profile model, and vitrinite reflectance of source rock and porosity of the reservoir were used as correction parameters, and the hydrocarbon migration and accumulation process was simulated by the invasion percolation algorithm. It is confirmed that the Ordovician Majiagou formation in the study area has three periods of hydrocarbon charging and three periods of hydrocarbon adjustment. The hydrocarbon migration pathways in carbonate rocks include lateral contact in the pinch-out area, source rock-reservoir contact in the trough area, fault-connected, and source rock-reservoir superposition. The source rock-reservoir-cap assemblages of upper-generation assemblage and lateral variable assemblage correspond to the lateral contact in the pinch-out area, source rock-reservoir contact in the trough area, fault-connected, while self-generation-reservoir-cap assemblage corresponds to the fault-connected, and source rock-reservoir superposition. The high part of the structure and fault zone are the strongest structural deformation parts, whose fractures are relatively developed, controlling the location of hydrocarbon accumulation, meanwhile, the reservoir quality directly controls the degree of hydrocarbons saturation. The research results are beneficial to establish a more accurate quantitative evaluation method of hydrocarbon migration and accumulation in carbonate rocks and to understanding the controlling factors of hydrocarbon migration and accumulation in carbonate rocks.

Effect of Thickness and Outlet Area Fraction of Macroporous Gas Diffusion Layers on Oxygen Transport Resistance in Water Injection Simulations

Enhanced water removal through the gas diffusion layer (GDL) is important for the design of high-performance proton exchange fuel cells. In this work, the effects of GDL thickness and open area fraction at the GDL/flow field interface are examined under water invasion for a carbon-paper GDL (similar to Toray TGP-H series). Both uncompressed and inhomogeneously compressed samples are considered. Transport in heterogeneous, macroporous GDLs is modeled by means of a hybrid 3D discrete/continuum formulation based on a subdivision of the porous medium into control volumes due to the lack of a well-defined separation between pore and layer scales. Capillary-dominated transport of liquid water is simulated with an invasion percolation algorithm, while oxygen diffusion is simulated with a continuum formulation. Model predictions are validated with previous numerical and experimental data. It is shown that the combination of thin GDLs (mathrm {thickness} sim 100; mu mathrm {m}) and high GDL/flow field open area fractions can facilitate water removal/oxygen supply from/to the catalyst layer and can provide a more uniform oxygen distribution over large cell active areas. In agreement with previous work, porous flow fields with pore sizes comparable to the GDL thickness are good candidates to meet the above requirements, while improving water removal from the flow field (higher gas-phase velocity than conventional millimeter-sized channels) and ensuring a more uniform assembly compression.

Pore network modeling of phase distribution and capillary force evolution during slow drying of particle aggregates

The spatial distribution of the liquid phase plays a crucial role for the damage and the deformation of particle aggregates during convective drying. In this paper, a three-dimensional triangulation pore network model (TPNM) is developed that can readily predict the global drying kinetics of particle aggregates and simulate the liquid phase distribution over time in the void space of aggregates made from spherical primary particles with different size distributions in the micrometer range. The discrete element method (DEM) is used to numerically generate the aggregates, whereby the complementary void space is decomposed into pores and throats using the regular Delaunay triangulation and its dual Voronoi tessellation. A modified version of the classical invasion percolation algorithm is set up to describe the preferential evaporation of the liquid confined in the pores. The liquid phase distribution obtained from isothermal drying TPNM simulations captures the evolution of phase patterns measured experimentally better than the prediction of the classical pore network model (CPNM) and the ring pore network model (RPNM). Local capillary forces caused by both fluid pressure and surface tension are computed over time from the filling state of pores. Taking a monodisperse particle aggregate as a reference, the size and spatial distributions of the primary particles are varied to represent the change in the corresponding void space morphology. The simulation results indicate the impact of the pore-scale heterogeneity on the spatial distribution of the liquid phase as well as the capillary forces during convective drying of particle aggregates.

A pore network approach to estimate the soil water retention curve of a granular material

Pore network modelling (PNM) is used to estimate the soil water retention curve (SWRC) of a granular soil, based on its grading properties. The ball and stick geometry is adopted, which considers the pores to be spherical, and the throats tubular. The topological characteristics of this geometry are evaluated based on the grain size distribution of the medium, and its grading properties. The pore size distribution is constructed based on the geometrical properties of the grading, and approached by a statistical distribution (normal or Weibull distributions). Additional hypotheses are proposed to establish the throat size distribution based on the neighboring pores. An invasion percolation algorithm is then applied to a saturated pore network, to simulate the introduction of air and the expulsion of water, thus mimicking the desorption of the sample. When analyzing the test data of uniformly graded sand, a reasonable agreement is observed between the experimental data and the numerical predictions.

A pore-network study on the factors influencing the isothermal drying of single- and dual-scale porous media

This paper analyzes a pore-network simulation to model isothermal drying in single-porosity and dual-porosity porous media after employing the invasion-percolation algorithm and the effect of liquid film. Utilizing a novel logistic equation to relate the inside and outside variables of the pore-network, the fully coupled implicit method is employed to solve the vapor transport problems in the pore-network and in the adjoining external flow field. A parametric study is conducted on the drying of a single-porosity square network placed next to a slit flow. Higher hydraulic diameter of the throats leads to higher drying rates and longer constant drying rate periods. The drying time increases and the drying rate decreases as the throat cross-section changes from a triangle to a square to a hexagon to a circle due to the weakening of the film effect. Increasing the external air humidity from 30% to 70% leads to a large decrease in the drying rates and the consequent increase in the overall drying times. Increasing the external flow velocity through Péclet number has little effect on the drying rates and times. Increasing the contact angle makes it harder for the corner films to materialize, which in turn decreases the constant drying rate period. Increase in the surface tension leads to an increase in the spread of the corner film which results in longer constant drying rate periods. Later the drying simulations of a dual-porosity domain using networks are carried out in a uniform air flow after keeping either the large pores or the small pores side open. The former leads to faster drying characterized by a previously observed bi-linear drying curve, and complete emptying of the outer large pores before the inner small pores. The latter witnesses the phenomenon of capillary pumping leading to a premature emptying of the inner larger pores. The case of all side open is also studied. Changing the throat cross-section from circle to square leads to a much faster drying due to the appearance of the film. Introduction of microstructural irregularity in the network by randomly changing the throat diameter and coordination number of pores does not affect the drying rate and drying time significantly.

Properties of the Vascular Networks in Malignant Tumors.

This work presents an analysis for real and synthetic angiogenic networks using a tomography image that obtains a portrait of a vascular network. After the image conversion into a binary format it is possible to measure various network properties, which includes the average path length, the clustering coefficient, the degree distribution and the fractal dimension. When comparing the observed properties with that produced by the Invasion Percolation algorithm (IPA), we observe that there exist differences between the properties obtained by the real and the synthetic networks produced by the IPA algorithm. Taking into account the former, a new algorithm which models the expansion of an angiogenic network through randomly heuristic rules is proposed. When comparing this new algorithm with the real networks it is observed that now both share some properties. Once creating synthetic networks, we prove the robustness of the network by subjecting the original angiogenic and the synthetic networks to the removal of the most connected nodes, and see to what extent the properties changed. Using this concept of robustness, in a very naive fashion it is possible to launch a hypothetical proposal for a therapeutic treatment based on the robustness of the network.

Pore network modeling of liquid water and oxygen transport through the porosity-graded bilayer gas diffusion layer of polymer electrolyte membrane fuel cells

In polymer electrolyte membrane fuel cells (PEMFCs), understanding the effect of graded porosity on liquid water distribution and oxygen diffusion inside the porous components is essential for the improvement of water management strategy. In this work, a regular three-dimensional pore network model is constructed to represent the bilayer gas diffusion layer (GDL) consisting of a porosity-graded micro-porous layer (GMPL) and a gas diffusion backing layer (GDBL). Based on this model, the liquid water distribution at breakthrough, relative oxygen effective diffusivity and limiting current density as a function of water saturation are obtained using invasion percolation algorithm and resistor-network theory, respectively. Parametric studies are also conducted to elucidate the effects of several structural factors, such as the graded porosity value of the GMPL, the GMPL thickness and the number of sub-layers in GMPL. The simulated results illustrate that the introduction of graded porosity can relieve the extent of the flooding phenomenon to a certain extent, and the underlying mechanism can be explained from the perspective of the pore scale. Furthermore, an optimum graded porosity value where the material has the preferred oxygen transfer capacity and high limit current density under partially saturated condition (water saturation is between 0.025 and 0.5) can be found.

Pore network modeling of fibrous porous media of uniform and gradient porosity

Two phase flow phenomena in fibrous media is studied using pore network modeling and visualized in 3D space. To explore the effect of network structure on flow property, virtual 3D fibrous mat with uniform and gradient porosity are constructed respectively using the realistic description of the complex fiber structure obtained from X-ray tomographic imaging. A marker-based watershed method is presented to extract pore networks directly from 3D segmented images. Invasion percolation algorithm is then applied to the extracted pore network to simulate primary drainage. The capillary pressure curve generated from drainage simulation is compared with an image-based full morphology method and achieves good agreement. The absolute and relative permeability are performed. For absolute permeability, we detected that the anisotropy of permeability tensor of layered fibrous mat can't be accurately reflected when pore cells are simplified using inscribed maximal balls. For relative permeability curves plotted as functions of water saturation, the results indicate a significant impact of the flow direction on liquid water distribution for porous material with gradient porosity.

Investigating the Effect of Non-Uniform Microporous Layer Intrusion on Oxygen Transport in Dry and Partially Saturated Polymer Electrolyte Membrane Fuel Cell Gas Diffusion Layers

State-of-the-art gas diffusion layers (GDL) used in polymer electrolyte membrane (PEM) fuel cells have a bi-layered structure comprising of a carbon fiber substrate and a microporous layer (MPL). Addition of the MPL has been shown to improve the fuel cell performance [1-2]. Earlier studies have also shown that the configuration of the MPL relative to the GDL has an effect on gas transport through the combined diffusion medium [3-4]. This study investigates the effect of the MPL intrusion and morphology on the oxygen transport in PEM fuel cell GDLs. In particular, we investigate the effect of varying MPL layer characteristics, both interior and exterior to the carbon fiber substrate in dry and saturated conditions, on the oxygen diffusivity and liquid water saturation at breakthrough conditions. Stochastic generation methods were utilized to create GDLs with varying MPL intrusion depths and morphologies. The generated morphologies of the non-uniform MPL layers were validated using reconstructed microscale X-ray-computed tomography. Pore network modelling was applied to simulate the oxygen transport behavior within these stochastically generated materials. The liquid water saturation was calculated using an invasion percolation algorithm. Three distinct scenarios were studied: varied MPL intrusion depth, varied MPL outer layer thickness, and constant MPL thickness with varied intrusion depth. In controlling the MPL intrusion, it was found that ~50-60 µm is the critical MPL thickness beyond which the saturated oxygen effective diffusion coefficient rapidly increases and the substrate saturation rapidly decreases (Figure 1). This study suggests that when the GDL is invaded by liquid water through cracks in the MPL, preferential MPL configurations exist for improving oxygen transport and minimizing substrate saturation. The results of this study can be used to guide the continued improvement of MPL fabrication.

Pore network modelling of molecular diffusion in a single‐block model during lean gas injection, investigating the effect of throat sorting

Based on analogy between isothermal drying of porous media and molecular diffusion in a single‐block model containing a volatile liquid during lean gas injection, a detailed procedure is presented for two‐dimensional quasi‐static pore network modelling of diffusion. A modified invasion percolation algorithm is used to account for evaporation of liquid clusters formed during gas injection process. Diffusion in gas phase, capillary‐induced flow in liquid phase, and evaporation in gas‐liquid interface were modelled to trace the movement of gas‐liquid interface and desaturation profile. Using a set of throat sizes with normal distribution, 15 cases of a regular pore network model were proposed by rearranging throat positions in the model. Results of this study indicate that diffusion time, breakthrough time, vapour partial pressure gradient, liquid saturation gradient, liquid cluster frequency, gas‐liquid interface frequency, and evaporation potential are strong functions of throat sorting. Correlation coefficient of throat size with position was used to capture the disorder of porous medium. It was shown that such a correlation coefficient is not capable of describing diffusion characteristics.

Quasi-static drainage in a network of nanoslits of non-uniform depth designed by grayscale laser lithography

A method is reported to fabricate silicon-glass nanofluidic chips with non-uniform channel depths in the range 20-500 nm and micrometer resolution in width. The process is based on grayscale laser lithography to structure photoresist in 2.5 dimensions in a single step, followed by a reactive ion etching to transfer the resist depth profile into silicon. It can be easily integrated in a complete process flow chart. The method is used to fabricate a network of interconnected slits of non-uniform depth, a geometry mimicking a nanoporous medium. The network is then used to perform a pressure step controlled drainage experiment, i.e. the immiscible displacement of a wetting fluid (liquid water) by a non-wetting one (nitrogen). The drainage patterns are analyzed by comparison with simulations based on the invasion percolation algorithm. The results indicate that slow drainage in the considered nanofluidic system well corresponds to the classical capillary fingering regime.

Achieving the Inside–Outside Coupling During Network Simulation of Isothermal Drying of a Porous Medium in a Turbulent Flow

The problem of isothermal (slow) drying of a porous medium in the form of a square bluff body resting on a flat surface while being exposed to turbulent air flow is investigated in this paper. The porous medium is represented by a pore-network model consisting of a network of volume-less nodes connected to each other via narrow throats. As the air flows past the water-saturated network, the water evaporates in the network and the evaporated water vapors diffuse through the air inside the network to advect away in the outside domain. In the history of such simulations coupling the ‘outside–inside’ processes during drying, this paper investigates, for the first time, the effect of turbulent flow on drying of a porous medium after coupling the outside-the-network flow and transport with the inside-the-network drying and liquid redistribution. Also, a commercial software package is used for the first time to obtain the outside flow. The water redistribution inside the network is predicted by the invasion-percolation algorithm after assuming the dominance of capillary forces over viscous and gravity forces. The outside velocity field is first obtained using the $$k\hbox {-}\varepsilon $$ model in ANSYS Fluent package. Then, the velocity field is used in a finite-volume-based code in order to solve for mass transfer outside the network using the regular convective-diffusive species transport equation. The drying mechanism and vapor transport inside the pore network is coupled with the outside mass-transfer simulation by considering a flux-balancing condition at the top surface of the network. After achieving grid independence and demonstrating the suitability of the no-slip boundary condition for the network top, the effects of the outside-flow Reynolds number and different turbulent flow models on several global drying parameters such as evaporation or drying rate, cumulative drying time, and top-surface saturation were studied. The resistance to vapor transport within the network was observed to dominate the enhanced vapor transport outside the network through the use of turbulent flow. Though the effect of variations in the network throat-size distribution upon saturation distribution within the network was found to be significant, its effect on the drying rate evolution was minimal. The concentration boundary layer thickness, which is employed to set the mass-transfer boundary condition in pore-network simulations and has been hitherto taken as constant in the literature, was found to change not only with position on the network top, but was also found to decline with time due to the drying of the network top.

Pore Network Modeling of Compressed Fuel Cell Components with OpenPNM

AbstractPore network modeling is used to model water invasion and multiphase transport through compressed PEFC gas diffusion layers. Networks are created using a Delaunay tessellation of randomly placed base‐points setting the pore locations and its compliment, the Voronoi diagram, is used to define the location of fibers and resultant pore and throat geometry. The model is validated in comparison to experimental capillary pressure curves obtained on compressed and uncompressed materials. Primary drainage is simulated with an invasion percolation algorithm that sequentially invades pores and throats separately with excellent agreement to experimental data, but required a slight modification to account for the higher aspect ratio of compressed pores. Compression is simulated by scaling the through‐plane coordinates in a uniform manner representing a GDL wholly beneath the current‐collector land. The relative permeability and diffusivity show some dependence on uniform compression. In‐plane porosity variations introduced by land‐channel compression are also investigated which have a marked effect on the limiting current. Saturation at breakthrough does not appear to be dependent on compression. However, a more important parameter, namely the peak saturation, is shown to influence the fuel cell performance and is dependent on the percolation inlet conditions.

Role of the microporous layer in the redistribution of phosphoric acid in high temperature PEM fuel cell gas diffusion electrodes

In this work, pore network modelling was used to investigate the role of the microporous layer (MPL) in the redistribution of phosphoric acid (PA) in a high temperature polymer electrolyte membrane fuel cell (HT-PEMFC) gas diffusion electrode (GDE). A HT-PEMFC GDE composed of a catalyst layer (sprayed in-house), MPL, and fibrous substrate was assembled, and its three-dimensional (3D) geometry was imaged using synchrotron X-ray microcomputed tomography. The spatial distribution of pore spaces and their connections were identified based on the 3D geometry, and an equivalent pore network of the GDE was obtained for simulating the PA transport with an invasion percolation algorithm. The predicted mass redistribution of PA in the GDE was found to be in excellent agreement with experimental values reported in the literature. The presence of the MPL was found to encourage the redistribution of PA within the CL and inhibited PA leaching toward the channel. To the authors’ best knowledge, this study is the first numerical simulation of PA transport in a HT-PEMFC GDE.

Impact of sample geometry on the measurement of pressure‐saturation curves: Experiments and simulations

AbstractIn this paper, we study the influence of sample geometry on the measurement of pressure‐saturation relationships, by analyzing the drainage of a two‐phase flow from a quasi‐2‐D random porous medium. The medium is transparent, which allows for the direct visualization of the invasion pattern during flow, and is initially saturated with a viscous liquid (a dyed glycerol‐water mix). As the pressure in the liquid is gradually reduced, air penetrates from an open inlet, displacing the liquid which leaves the system from an outlet on the opposite side. Pressure measurements and images of the flow are recorded and the pressure‐saturation relationship is computed. We show that this relationship depends on the system size and aspect ratio. The effects of the system's boundaries on this relationship are measured experimentally and compared with simulations produced using an invasion percolation algorithm. The pressure build up at the beginning and end of the invasion process are particularly affected by the boundaries of the system whereas at the central part of the model (when the air front progresses far from these boundaries), the invasion happens at a statistically constant capillary pressure. These observations have led us to propose a much simplified pressure‐saturation relationship, valid for systems that are large enough such that the invasion is not influenced by boundary effects. The properties of this relationship depend on the capillary pressure thresholds distribution, sample dimensions, and average pore connectivity and its applications may be of particular interest for simulations of two‐phase flow in large porous media.

RELATIVE PERMEABILITIES CHARACTERIZATION IN CHEMICAL FLOODING WITH THE CONSIDERATION OF VISCOSITY RATIO AND INTERFACIAL TENSION BY A PORE-SCALE NETWORK MODEL

A dynamic pore-scale network model is presented for investigating the effects of interfacial tension and oil−water viscosity on relative permeability during chemical flooding. This model takes into account both viscous and capillary forces in analyzing the impacts of chemical properties on flow behavior or displacement configuration, as opposed to the conventional or invasion percolation algorithm which incorporates capillary pressure only. The study results indicate that both water and oil relative-permeability curves are dependent strongly on interfacial tension as well as an oil−water viscosity ratio. In particular, water and oil relative-permeability curves are both found to shift upward as interfacial tension is reduced, and they both tend to become linear versus saturation once interfacial tension is at low values. In addition, the oil−water viscosity ratio appears to have only a small effect under conditions of high interfacial tension. When the interfacial tension is low, however, water relative permeability decreases more rapidly (with the increase in the aqueous-phase viscosity) than oil relative permeability. The breakthrough saturation of the aqueous phase during chemical flooding tends to decrease with the reduction of interfacial tension and may also be affected by the oil−water viscosity ratio.

Pore-network modeling of liquid water flow in gas diffusion layers of proton exchange membrane fuel cells

Abstract Fluid flow through the gas diffusion layer (GDL) of fuel cells is numerically studied using a pore network modeling approach. The model is developed based on an experimental visualization technique (fluorescence microscopy). The images obtained from this technique are analyzed to find patterns of flow inside the GDL samples with different hydrophobicity. Three different flow patterns are observed: initial invasion, progression, and pore-filling. The observation shows that liquid water flows into the majority of available pores on the boundary of the untreated GDL and several branches are segregated from the initial pathways. For the treated GDL, however, a handful of boundary pores are invaded and the original pathways extend toward the other side of the medium with minimum branching. The numerical model, developed based on an invasion percolation algorithm, is used to study the effects of GDL hydrophobicity and thickness on the flow configuration and breakthrough time as well as to determine the flow rate and saturation in different GDL samples. During the injection of water into the samples, it is numerically shown that the flow rates are monotonically decreasing for both treated and untreated samples. For the treated sample, however, the injection flow rate is constantly lower than that of the untreated sample, resulting in a lower overall water saturation at breakthrough. The numerical results also suggest that hydrophobic treatment of thick samples has minor effects on water management and overall performance. The developed model can be used to optimize the GDL properties for designing porous medium with effective transport characteristics.