- About this Journal
- Abstracting and Indexing
- Aims and Scope
- Annual Issues
- Article Processing Charges
- Articles in Press
- Author Guidelines
- Bibliographic Information
- Citations to this Journal
- Contact Information
- Editorial Board
- Editorial Workflow
- Free eTOC Alerts
- Publication Ethics
- Reviewers Acknowledgment
- Submit a Manuscript
- Subscription Information
- Table of Contents
Advances in Mechanical Engineering
Volume 2010 (2010), Article ID 142879, 11 pages
Mesoscopic Modeling of Multiphysicochemical Transport Phenomena in Porous Media
Computational Earth Science Group, Earth and Environmental Sciences Division, Los Alamos National Laboratory, Los Alamos, NM 87545, USA
Received 6 October 2009; Accepted 11 December 2009
Academic Editor: Chen Li
Copyright © 2010 Qinjun Kang et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
We present our recent progress on mesoscopic modeling of multiphysicochemical transport phenomena in porous media based on the lattice Boltzmann method. Simulation examples include injection of -saturated brine into a limestone rock, two-phase behavior and flooding phenomena in polymer electrolyte fuel cells, and electroosmosis in homogeneously charged porous media. It is shown that the lattice Boltzmann method can account for multiple, coupled physicochemical processes in these systems and can shed some light on the underlying physics occurring at the fundamental scale. Therefore, it can be a potential powerful numerical tool to analyze multiphysicochemical processes in various energy, earth, and environmental systems.
Multiphysicochemical transport phenomena in porous media are ubiquitous, particularly in various energy, earth, and environment systems. One example is the disposal of supercritical CO2 in geologic formations, the most promising near-term solution to the problem of reducing carbon emissions into the atmosphere . Experimental analyses of the long-term fate of CO2 after injection into the geologic formations are not possible with relatively short-term laboratory experiments. Therefore it is necessary to employ comprehensive numerical models that incorporate multiple physicochemical processes involving interactions between the injected CO2, the brine in the pore spaces, and the minerals lining the pores. Supercritical CO2, as a buoyant and slightly miscible fluid, once injected, displaces brine from the pore space in a complex pattern. At the interface with brine, CO2 dissolves into the brine to form carbonic acid that can react with and dissolve minerals eventually leading to mineral precipitation further along the flow path. Clearly, there are multiple physics processes involved, including hydrodynamics, thermodynamics, chemical dynamics, and electrodynamics (because the surface of most natural media is charged). All these processes are ultimately governed by pore-scale interfacial phenomena, which occur at scales of microns. However, because of the wide disparity in scales ranging from pore to field, a continuum formulation based on spatial averages taken over length scales much larger than typical pore and mineral grain sizes is often utilized, implying the existence of a representative elemental volume (REV) . As a result, spatial heterogeneities at smaller scales are unresolved and the aggregate effects of the pore-scale (mesoscopic scale) processes are taken into account through various effective constitutive parameters. One of the goals of performing pore-scale simulations is to obtain values for these constitutive parameters through upscaling the pore-scale results. Other goals are to identify key parameters and physicochemical processes that control macroscopic phenomena, and to validate continuum descriptions.
Another example is fuel cells, and in particular polymer electrolyte fuel cells (PEFCs). In PEFCs, the catalyst layer (CL) is the host to several competing transport mechanisms involving charge (proton and electron), species (oxygen, nitrogen, and water vapor), and liquid water transport. The multi-faceted functionality of a gas diffusion layer (GDL) includes reactant distribution, liquid water transport, electron transport, heat conduction, and mechanical support to the membrane-electrode-assembly. Despite tremendous recent progress in enhancing the overall cell performance, a pivotal performance limitation in PEFCs is manifested in terms of mass transport loss originating from suboptimal liquid water transport and resulting flooding in the constituent components . In recent years, several macroscopic computational models for multiple-physicochemical transport processes in PEFCs [4–10] have been developed. These macroscopic models again are based on the theory of volume averaging and treat the catalyst layer and gas diffusion layer as macrohomogeneous porous layers. Due to their macroscopic nature, the current models fail to resolve the influence of the structural morphology of the CL and GDL on the underlying physics. Mesoscopic modeling is critical to understanding of the overall structure-wettability-transport interactions as well as the underlying multiphysicochemical processes in the CL and GDL, and hence is a useful tool for design and optimization of microstructures of electrodes for better performance and durability.
In this paper, we review our recent work on mesoscopic modeling of multiphysicochemical processes in porous media, based on the lattice Boltzmann method (LBM), a relatively new numerical method for simulating fluid flows and modeling physics in fluids . Originating from the classical statistical physics, LBM is a mesoscopic method based on simplified kinetic equations. In the LBM, the fluid is modeled as a collection of fictitious particles propagating and colliding over a discrete lattice domain. Mesoscopic continuity and momentum equations can be obtained from this propagation-collision dynamics through a rigors mathematical analysis. The particulate nature and local dynamics provide advantages for complex boundaries and parallel computation. In addition, the kinetic nature of the LBM makes it easy to account for new physics in the LBM framework, which is especially useful for modeling multiphysicochemical phenomena. In Section 2, the partial differential equations governing fluid flow, transport of reactive species and electric potential, as well as mineral reactions in porous media will be given. In Section 3, the implementation of the LBM to solve these governing equations will be presented. Some simulation examples will be given in Section 4 and concluding remarks in Section 5.
2. Governing Equations
Consider the electrolyte fluid flowing through solid porous media. Although the pore scale in this study may be at micrometer scale, the fluid can still be treated as a continuum Newtonian fluid since the characteristic size is much larger than the molecular diameters [12–14].
2.1. Continuity and Momentum Equations
For isothermal incompressible fluid flow, the continuity and momentum equations can be written as 
where represents the density of the fluid, the time, the velocity vector, the pressure, the kinetic viscosity, and the body force density which may include all the effective body forces.
For the hydrodynamic boundary condition, we use the nonslip model at the solid surfaces. The slip boundary conditions have been adopted in some recent studies [16, 17]; however, a careful molecular study showed that the hydrodynamic boundary condition, slip or not, depended on the molecular interactions between fluid and solid and on the channel size [18, 19]. For the flow in porous media considered in this work, the nonslip boundary condition is still valid.
2.2. Transport Equations for Aqueous Species and Electrical Potential
For the th ion species in the solute, the mass conservation equation describing transport and reaction can be written in the general form 
where denotes the ionic concentration, the species flux, a radioactive decay constant, and the rate at which the ith species is produced or consumed by chemical reactions. The flux , consisting of contributions from advection, diffusion, and electrochemical migration terms has the form 
where the first term on the right refers to electrochemical migration, the second term to aqueous diffusion, and the last term to advective transport. Here , , and denote the ion algebraic valence, the diffusivity, and the activity coefficient of the th species, respectively; and , , and denote the absolute charge of electron, the Boltzmann constant, and the absolute temperature, respectively. The quantity represents the local electrical potential caused by the ionic distribution which is governed by the Poisson equation
where is the local dimensionless fluid dielectric constant, the permittivity of a vacuum, and the net charge density. Assuming no radiation and constant activity coefficient and substituting (4) into (2), we have
This is the Nernst-Planck equation , where can be any kind of effective body force. In this contribution we only consider the static electrical force from an external electric field. The general form of electrical force in electrokinetic fluids can be expressed as
where is the external electrical potential field.
When the ionic convection is negligible and the electric potential is continuously derivable, (5) has a simple steady-state solution for dilute electrolyte solutions:
2.3. Equations for Mineral Reaction Rates
Heterogeneous reactions between aqueous species and minerals at the pore-mineral interface are given by 
where denotes the unit normal perpendicular to the fluid-mineral interface pointing toward the fluid phase, D denotes the aqueous diffusion coefficient assumed to be the same for all species for simplicity, and denotes the reaction flux of the th mineral at its surface, taken as positive for precipitation and negative for dissolution. The total surface area across which diffusion takes place equal to the sum of individual mineral surface areas is given as
3. Lattice Boltzmann Model Implementation
3.1. Incompressible Lattice Boltzmann Model for Single-Phase Flow
In order to eliminate compressible effects in the conventional LBM, here we use an incompressible LB model constructed by Guo et al. . The pore-scale flow of a single aqueous fluid phase is simulated by the following evolution equation:
In the above equation, is the time increment, the distribution function along the direction in velocity space, the corresponding equilibrium distribution function, and τ the dimensionless relaxation time. For the commonly used two-dimensional, nine-speed LB model (D2Q9), the discrete velocities have the following form:
For the incompressible LB model, the equilibrium distribution is defined by Guo et al. :
where , and are the parameters satisfying
In these equations, , where is the space increment, and and are the pressure and velocity of the fluid, respectively. The corresponding weight coefficients are
3.2. Lattice Boltzmann Model for Two-Phase Flow
The interaction-potential model, originally proposed by Shan and Chen , and henceforth referred to as the S-C model, introduces distribution functions for a fluid mixture comprising of components. Each distribution function represents a fluid component and satisfies the evolution equation. The nonlocal interaction between particles at neighboring lattice sites is included in the kinetics through a set of potentials. The LB equation for the th component can be written as
where is the number-density distribution function for the th component in the th velocity direction at position and time , and is the time increment. In the term on the right-hand side, is the relaxation time of the kth component in lattice unit, and is the corresponding equilibrium distribution function.
The phase separation between different fluid phases, the wettability of a particular fluid phase to the solid, and the body force are taken into account by modifying the velocity used to calculate the equilibrium distribution function. An extra component-specific velocity due to interparticle interaction is added on top of a common velocity for each component. Interparticle interaction is realized through the total force, , acting on the kth component, including fluid/fluid interaction, fluid/solid interaction, and external force. More details can be found in [26, 27].
The continuity and momentum equations can be obtained for the fluid mixture as a single fluid using Chapman-Enskog expansion procedure in the nearly incompressible limit:
where the total density and velocity of the fluid mixture are given, respectively, by
with a nonideal gas equation of state given by .
3.3. Lattice Boltzmann Model for Transport of Reactive Solutes
In a previous article, Kang et al.  have derived the following LB equation for the total primary species concentrations for chemical systems with reactions written in canonical form:
where is the number of primary species, is the total concentration of the th primary species, is its distribution function along the α direction, is the corresponding equilibrium distribution function, are velocity vectors, is the time increment, and is the dimensionless relaxation time.
It has been shown that the above equation can recover the following pore-scale advection-diffusion equation for :
This equation is the same as (5) except that here the electrochemical migration is neglected. Assuming that the homogeneous reactions are in instantaneous equilibrium, we have the following mass action equation [31, 32]:
where are the stoichiometric coefficients, is the equilibrium constant of the ith homogeneous reaction, is the activity coefficient of the th secondary species, and and are solute concentrations for primary and secondary species, respectively. They are related by
where is the number of independent homogeneous reactions, or, equivalently, secondary species.
3.4. Lattice Poisson Method
To solve the Poisson equation with strong nonlinearity, (8), we adopt the lattice Poisson method (LPM) developed previously [34, 35], which tracks the electrical potential distribution transporting on the discrete lattices. By expanding (8) into the time-dependent form
with representing the negative right-hand side (RHS) term of the original (8), we get the discrete evolution equation for the electrical potential distribution
where is the equilibrium distribution of the electric potential evolution variable. The time step for the electrical potential evolution is
where is a pseudo-sound speed in the potential field. After evolving on the discrete lattices, the mesoscopic electrical potential can be calculated using
Although the electrical potential evolution equations are in an unsteady form, only the steady state result is realistic, because the electromagnetic susceptibility has not been considered. Although the lattice evolution method for nonlinear Poisson equation is not as efficient as the multi-grid solutions due to its long wavelength limit, it has the advantages of suitability for geometrical complexity and parallel computing.
4. Simulation Examples
4.1. Injection of CO2 into a Limestone Rock
We first present some modeling results on the injection of a fluid saturated with 170 bars CO2(g) into a limestone rock at the pore scale. The pore structure was derived from a digitized image of a limestone rock thin section with pixels (Figure 1). We reduced the original resolution to save computational time (Figure 2). The chemical system of Na+-Ca2+-Mg2+-H+--Cl--CO2 with the reaction of calcite to form dolomite and gypsum is considered. Secondary species included in the simulation are OH-, , H2SO4(aq), , , CaCO3(aq), , CaOH+, CaSO4(aq), MgCO3(aq), , MgSO4(aq), MgOH+, NaCl(aq), NaHCO3(aq), and NaOH(aq). Initial fluid composition is pH 7.75 and 2.69 m NaCl brine, equilibrium with minerals calcite, dolomite and gypsum at 25°C. Initial rock composition is calcite. Secondary minerals include dolomite and gypsum. For boundary conditions, a constant pressure gradient is imposed across the domain for flow. When flow reaches steady state, a fluid with a pH of 3.87 and in equilibrium with 179 bars CO2(g) and minerals dolomite and gypsum is introduced at the inlet. Zero gradient boundary conditions are imposed at the outlet. Two different cases are considered with different mineral reaction rates to show their effects on solution concentration, mineral deposition, and change in geometry.
Resulting geometries at time = 15625 seconds for two different mineral reaction rate constants are plotted in Figure 3. Damkohler is 7.375 for calcite and gypsum, 0.7375 for dolomite for the faster mineral reactions, for calcite and gypsum, and for dolomite for slower reactions. Concentration distribution of total Ca2+, Mg2+, and , pH, and reaction rates for calcite, dolomite, and gypsum for the slower reactions are plotted in Figure 4. As can been seen from the figures, as the reaction rate constants decrease, the deposition of dolomite becomes more uniform surrounding the dissolving calcite grains. Only a small amount of gypsum forms on top of dolomite. At some point in the simulation, the major pores for flow become blocked halting further fluid flow through the medium. The pH is uniform over the entire pore space. All reaction rates have finite values at the mineral surface in the whole domain, outlining the solid geometry. The reaction rate is negative for calcite and positive for dolomite and gypsum, confirming that calcite is dissolving while dolomite and gypsum are precipitating.
4.2. Two-Phase Behavior and Flooding Phenomena in Polymer Electrolyte Fuel Cells
In this Section, we present some results for two-phase flow through the porous CL and the fibrous GDL in a PEFC. Details can be found in . Figure 5 displays the steady state invading liquid water fronts corresponding to increasing capillary pressures from the primary drainage simulation in the reconstructed CL microstructure characterized by slightly hydrophobic wetting characteristics with a static contact angle of . At lower capillary pressures, the liquid water saturation front exhibits finger like pattern, similar to the displacement pattern observed typically in the capillary fingering regime. The displacing liquid water phase penetrates into the body of the resident wetting phase (i.e., air) in the shape of fingers owing to the surface tension driven capillary force. However, at high saturation levels, the invading nonwetting phase tends to exhibit a somewhat flat advancing front. This observation, as highlighted in Figure 5(b), indicates that with increasing capillary pressure, even at very low capillary number (Ca), several penetrating saturation fronts tend to merge and form a stable front. The invasion pattern transitions from the capillary fingering regime to the stable displacement regime and potentially lies in the transition zone in between. In an operating fuel cell, the resulting liquid water displacement pattern pertaining to the underlying pore-morphology and wetting characteristics would play a vital role in the transport of the liquid water and hence the overall flooding behavior.
Figure 6 shows the liquid water distribution as well as the invasion pattern from the primary drainage simulation with increasing capillary pressure in the initially air-saturated reconstructed carbon paper GDL characterized by hydrophobic wetting characteristics with a static contact angle of . The reconstructed GDL structure used in the two-phase simulation consists of lattice points in order to manage the computational overhead to a reasonable level. Physically, this scenario corresponds to the transport of liquid water generated from the electrochemical reaction in a hydrophobic CL into the otherwise air-occupied GDL in an operating fuel cell. At the initially very low capillary pressure, the invading front overcomes the barrier pressure only at some preferential locations depending upon the pore size along with the emergence of droplets owing to strong hydrophobicity. As the capillary pressure increases, several liquid water fronts start to penetrate into the air occupied domain. Further increase in capillary pressure exhibits growth of droplets at two invasion fronts, followed by the coalescence of the drops and collapsing into a single front. This newly formed front then invades into the less tortuous in-plane direction. Additionally, emergence of tiny droplets and subsequent growth can be observed in the constricted pores in the vicinity of the inlet region primarily due to strong wall adhesion forces from interactions with highly hydrophobic fibers with the increasing capillary pressure. One of the several invading fronts finally reaches the air reservoir, physically the GDL/channel interface, at a preferential location corresponding to the capillary pressure and is also referred to as the bubble point. It is important to note that the mesoscopic LB simulations provide fundamental insight into the pore-scale liquid water transport through different GDL structures and would likely enable novel GDL design with better water removal and flooding mitigation.
4.3. Electroosmosis in Homogeneously Charged Micro- and Nanoscale Porous Media
In this section, we briefly present some simulation results on electroosmotic flows (EOFs) in charged micro porous media using the lattice Poisson-Boltzmann method (LPBM), with geometry effects, solution and surface charge effects considered. More details can be found in . We focus on a cubic system of which each side is 1 micron long. A uniform grid is used. We change microstructure geometries of porous media by varying the porosity from 0.1 to 0.9. The average characteristic length of particles varies from 20 to 150 nm. Figure 7 shows Schematics of the generated porous structures for porosity 0.3 and 0.6. The bulk ionic concentration varies from 10-6 to 10-3 M and the surface zeta potential from 0 to mV. The other properties and parameters used in this work are the fluid density kg/m3, the dielectric constant C2/Jm, the dynamic viscosity mPas, the temperature K and the external electrical field strength V/m.
First, the geometry effects on the electroosmotic permeability in micro porous media are investigated by changing volume fraction and particle size (or number density) of the solid phase. We define the electroosmotic permeability, , as
where is the averaged velocity of EOF along the direction of the driving electrical field . The coefficients of electroosmotic permeability () for different porosities () of porous media are shown in Figure 8. The bulk molar concentration M, and mV. The electroosmotic permeability increases with the porosity monotonically. The increasing rate rises with the porosity as well which is very low when the porosity is smaller than 0.5 and becomes sharply high when the porosity is larger than 0.7. The predicted electroosmotic permeability is in the order of 10-9 m2/sV, which is consistent with the existing experimental measurements.
Figure 9 shows the predicted electroosmotic permeability versus the bulk ionic concentration of the electrolyte solution. We used . The electroosmotic permeability increase monotonically with the bulk ionic concentration as varies from 10-6 to 10-3 M. This result can be explained by the undeveloped electrical potential distributions in narrow channels, whose similar results can be found in [37, Figure ] and [38, Figures and ]. When is lower than 10-4 M, the electroosmotic permeability is nearly proportional to the bulk ionic concentration. When is higher, the increasing rate becomes a little smaller.
Zeta potential on solid surfaces of porous media affects EOF permeability directly. Simple proportional relationships have been obtained between the electroosmotic permeability and the zeta potential for electrical transports in soils and in polymer composites recently based on the boundary-layer theory. Here we analyze such effects using our numerical methods.
Figure 10 shows the calculated electroosmotic permeability versus the zeta potential on solid surfaces of porous media. All surfaces are homogeneously charged with a same value of . The other parameters used are: M, and . The zeta potential changes from 0 to 90 mV. It shows that the proportionally linear relationship between electroosmotic permeability and zeta potential is accurate only when is very small (30 mV). The permeability increases much sharper when the zeta potential is larger than 40 mV.
We have presented our recent work on mesoscopic modeling of multiphysicochemical processes in porous media, based on the LBM. For the problem of injecting CO2 saturated brine into a limestone rock, it is shown that the LBM is able to provide detailed information on fluid velocity, solute concentration, mineral composition, and reaction rates, as well as the evolution of the porous media geometry, and therefore can shed some light on the fundamental physics occurring at the pore scale for reactive transport involved in geologic CO2 sequestration. For two-phase behavior and flooding phenomena in PEFCs, the LBM is a powerful tool to study the influence of the pore structure and surface wettability on liquid water transport and interfacial dynamics in the PEFC catalyst layer and gas diffusion layer. Particularly, the two-phase regime transition phenomenon in the capillary dominated transport in the CL and the influence of the mixed wetting characteristics on the flooding dynamics in the GDL are demonstrated. For electroosmotic flows in charged porous media, the strongly nonlinear governing equations of electroosmosis in three-dimensional porous media are solved using the LPBM. The effects of pore structure, bulk ionic concentration, and the surface charge on electroosmotic permeability are carefully investigated. It is concluded that the LBM is a powerful numerical tool to simulate multiphysicochemical processes in porous media at the pore-scale.
This article is based on research project supported by the National Science Foundation under Grant no. CHE-0431328 and the U.S. Department of Energy, Biological and Environmental Research (BER), and by LDRD projects 20070760PRD4 and 20080727PRD2 sponsored by Los Alamos National Laboratory.
- S. Pacala and R. Socolow, “Stabilization wedges: solving the climate problem for the next 50 years with current technologies,” Science, vol. 305, no. 5686, pp. 968–972, 2004.
- J. Bear, Dynamics of Fluids in Porous Media, Elsevier, New York, NY, USA, 1972.
- C.-Y. Wang, “Fundamental models for fuel cell engineering,” Chemical Reviews, vol. 104, no. 10, pp. 4727–4766, 2004.
- J. J. Baschuk and X. Li, “Modelling of polymer electrolyte membrane fuel cells with variable degrees of water flooding,” Journal of Power Sources, vol. 86, no. 1, pp. 181–196, 2000.
- L. You and H. Liu, “A two-phase flow and transport model for the cathode of PEM fuel cells,” International Journal of Heat and Mass Transfer, vol. 45, no. 11, pp. 2277–2287, 2002.
- W. He, J. S. Yi, and T. Van Nguyen, “Two-phase flow model of the cathode of PEM fuel cells using interdigitated flow fields,” AIChE Journal, vol. 46, no. 10, pp. 2053–2064, 2000.
- J. H. Nam and M. Kaviany, “Effective diffusivity and water-saturation distribution in single- and two-layer PEMFC diffusion medium,” International Journal of Heat and Mass Transfer, vol. 46, no. 24, pp. 4595–4611, 2003.
- S. Dutta, S. Shimpalee, and J. W. Van Zee, “Numerical prediction of mass-exchange between cathode and anode channels in a PEM fuel cell,” International Journal of Heat and Mass Transfer, vol. 44, no. 11, pp. 2029–2042, 2001.
- M. Noponen, E. Birgersson, J. Ihonen, M. Vynnycky, A. Lundblad, and G. Lindbergh, “A two-phase non-isothermal PEFC model: theory and validation,” Fuel Cells, vol. 4, no. 4, pp. 365–377, 2004.
- Z. Wang, X. Wu, R. Ni, and Y. Wang, “Binocular fusion in Panum's limiting case of stereopsis obeys the uniqueness constraint,” Science in China, Series C, vol. 44, no. 1, pp. 40–48, 2001.
- S. Chen and G. D. Doolen, “Lattice boltzmann method for fluid flows,” Annual Review of Fluid Mechanics, vol. 30, pp. 329–364, 1998.
- M. Wang and S. Chen, “On applicability of Poisson-Boltzmann equation for micro- and nanoscale electroosmotic flows,” Communications in Computational Physics, vol. 3, no. 5, pp. 1087–1099, 2008.
- M. Wang, J. Liu, and S. Chen, “Similarity of electroosmotic flows in nanochannels,” Molecular Simulation, vol. 33, no. 3, pp. 239–244, 2007.
- M. Wang, J. Liu, and S. Chen, “Electric potential distribution in nanoscale electroosmosis: from molecules to continuum,” Molecular Simulation, vol. 33, no. 15, pp. 1273–1277, 2007.
- L. D. Landau and E. M. Lifshitz, Fluid Mechanics, Elsevier, New York, NY, USA, 1959.
- C. Davidson and X. Xuan, “Electrokinetic energy conversion in slip nanochannels,” Journal of Power Sources, vol. 179, no. 1, pp. 297–300, 2008.
- Y. Ren and D. Stein, “Slip-enhanced electrokinetic energy conversion in nanofluidic channels,” Nanotechnology, vol. 19, no. 19, Article ID 195707, 2008.
- J.-F. Dufrêche, V. Marry, N. Malíková, and P. Turq, “Molecular hydrodynamics for electro-osmosis in clays: from Kubo to Smoluchowski,” Journal of Molecular Liquids, vol. 118, no. 1–3, pp. 145–153, 2005.
- L. Joly, C. Ybert, E. Trizac, and L. Bocquet, “Hydrodynamics within the electric double layer on slipping surfaces,” Physical Review Letters, vol. 93, no. 25, Article ID 257805, 4 pages, 2004.
- P. C. Lichtner, “Principles and practice of reactive transport modeling,” Materials Research Society Symposium Proceedings, vol. 353, no. 1, pp. 117–130, 1995.
- V. G. Levich, Physico-Chemical Hydrodynamics, Prentice-Hall, New York, NY, USA, 1962.
- B. Honig and A. Nicholls, “Classical electrostatics in biology and chemistry,” Science, vol. 268, no. 5214, pp. 1144–1149, 1995.
- P. C. Lichtner and Q. Kang, “Upscaling pore-scale reactive transport equations using a multiscale continuum formulation,” Water Resources Research, vol. 43, no. 12, Article ID W12S15, 2007.
- Z. Guo, B. Shi, and N. Wang, “Lattice BGK model for incompressible Navier-Stokes equation,” Journal of Computational Physics, vol. 165, no. 1, pp. 288–306, 2000.
- X. Shan and H. Chen, “Lattice Boltzmann model for simulating flows with multiple phases and components,” Physical Review E, vol. 47, no. 3, pp. 1815–1819, 1993.
- Q. Kang, D. Zhang, and S. Chen, “Displacement of a two-dimensional immiscible droplet in a channel,” Physics of Fluids, vol. 14, no. 9, pp. 3203–3214, 2002.
- Q. Kang, D. Zhang, and S. Chen, “Displacement of a three-dimensional immiscible droplet in a duct,” Journal of Fluid Mechanics, vol. 545, pp. 41–66, 2005.
- X. Shan and G. Doolen, “Diffusion in a multicomponent lattice Boltzmann equation model,” Physical Review E, vol. 54, no. 4, pp. 3614–3620, 1996.
- Q. Kang, P. C. Lichtner, and D. Zhang, “Lattice Boltzmann pore-scale model for multicomponent reactive transport in porous media,” Journal of Geophysical Research B, vol. 111, no. 5, Article ID B05203, 2006.
- S. P. Dawson, S. Chen, and G. D. Doolen, “Lattice Boltzmann computations for reaction-diffusion equations,” The Journal of Chemical Physics, vol. 98, no. 2, pp. 1514–1523, 1993.
- P. C. Lichtner, C. I. Steefel, and E. H. Oelkers, Eds., Reactive Transport in Porous Media, vol. 34 of Reviews in Mineralogy, P. H. Ribbe, Mineralogical Society of America, Washington, DC, USA, 1996.
- P. C. Lichtner, “Continuum model for simultaneous chemical reactions and mass transport in hydrothermal systems,” Geochimica et Cosmochimica Acta, vol. 49, no. 3, pp. 779–800, 1985.
- Q. Kang, P. C. Lichtner, and D. Zhang, “An improved lattice Boltzmann model for multicomponent reactive transport in porous media at the pore scale,” Water Resources Research, vol. 43, no. 12, Article ID W12S14, 2007.
- J. Wang, M. Wang, and Z. Li, “Lattice Poisson-Boltzmann simulations of electro-osmotic flows in microchannels,” Journal of Colloid and Interface Science, vol. 296, no. 2, pp. 729–736, 2006.
- M. Wang and S. Chen, “Electroosmosis in homogeneously charged micro- and nanoscale random porous media,” Journal of Colloid and Interface Science, vol. 314, no. 1, pp. 264–273, 2007.
- P. P. Mukherjee, C.-Y. Wang, and Q. Kang, “Mesoscopic modeling of two-phase behavior and flooding phenomena in polymer electrolyte fuel cells,” Electrochimica Acta, vol. 54, no. 27, pp. 6861–6875, 2009.
- M. Wang, N. Pan, J. Wang, and S. Chen, “Lattice Poisson-Boltzmann simulations of electroosmotic flows in charged anisotropic porous media,” Communications in Computational Physics, vol. 2, no. 6, pp. 1055–1070, 2007.
- M. Wang, J. Wang, and S. Chen, “Roughness and cavitations effects on electro-osmotic flows in rough microchannels using the lattice Poisson-Boltzmann methods,” Journal of Computational Physics, vol. 226, no. 1, pp. 836–851, 2007.