PHYSICAL REVIEW LETTERS PRL 101 044502 2008 week ending 25 JULY 2008 Stochastic Langevin Model for Flow and Transport in Porous Media Alexandre M Tartakovsky Pacific Northwest National Laboratory Richland Washington 99352 USA Daniel M Tartakovsky University of California San Diego La Jolla California 92093 USA Paul Meakin Idaho National Laboratory Idaho Falls Idaho 83415 USA Received 26 November 2007 revised manuscript received 27 April 2008 published 22 July 2008 We present a new model for fluid flow and solute transport in porous media which employs smoothed particle hydrodynamics to solve a Langevin equation for flow and dispersion in porous media This allows for effective separation of the advective and diffusive mixing mechanisms which is absent in the classical dispersion theory that lumps both types of mixing into dispersion coefficient The classical dispersion theory overestimates both mixing induced effective reaction rates and the effective fractal dimension of the mixing fronts associated with miscible fluid Rayleigh Taylor instabilities We demonstrate that the stochastic Langevin equation model overcomes these deficiencies DOI 10 1103 PhysRevLett 101 044502 PACS numbers 47 56 r 05 10 Gg Flow and transport in porous media can be described on two fundamental scales On the pore scale these phenomena are governed by the Navier Stokes and advectiondiffusion equations Detailed knowledge of the pore geometry of most natural and manufactured porous media is elusive and solving these equations over a large volume of a porous medium is often computationally prohibitive These and other practical considerations have led to the development of continuum or Darcy scale models which are obtained by averaging the Navier Stokes and or advection diffusion equations over a sufficiently large volume of the porous medium After a series of simplifying assumptions the volumetric or statistical averaging of the pore scale continuity and Navier Stokes equations yields 1 the Darcy scale continuity equation d r u dt 1 and the Darcy scale momentum conservation equation rp du g u dt 2 Here d dt t u r denotes the material derivative u is the mean microscopic velocity of a fluid with density and viscosity g is the gravitational acceleration and p is the pressure The friction coefficient k porosity permeability k and hydraulic conductivity g are some of the macroscopic parameters characterizing porous media on the Darcy scale While it is common to simplify 2 further by setting du dt to 0 which gives rise to Darcy s law we retain this term for completeness A similar averaging procedure applied to the pore scale advection diffusion equation gives rise to the Darcy scale 0031 9007 08 101 4 044502 4 advection dispersion equation ADE 1 dC r D rC dt 3 Here C is the solute concentration defined as a mass of solute per unit mass of a solution D Dm juj is the dispersion coefficient written as a scalar rather than a tensor to simplify the presentation Dm is the molecular diffusion coefficient and and are the tortuosity and the dispersivity of the porous medium While adequate in many settings the classical Darcyscale equations 1 3 have a number of known conceptual and operational drawbacks For example 3 implies that hydrodynamic dispersion is functionally analogous to Fickian diffusion with a macro scale effective diffusion coefficient D that lumps together advective mixing spreading due to variations in the fluid velocity and diffusive mixing 1 This contradicts a number of observations which indicate that dispersive mixing is fundamentally different from its purely diffusive counterpart Specifically the fractal dimensions of the diffusion and dispersion fronts isoconcentration contours are different 2 and ADE based models of reactive transport can significantly over predict the extent of reaction in mixinginduced chemical transformations 3 5 In heterogeneous porous media many of the shortcomings of the traditional ADE can be overcome by treating the hydraulic conductivity or equivalently the mean microscopic velocity u in 3 as random fields and applying either stochastic averaging e g 6 or renormalizationgroup analyses e g 7 According to these and other similar approaches the randomness is absent if the porous medium is homogeneous and the drawbacks of the ADE 3 reemerge 044502 1 2008 The American Physical Society PRL 101 044502 2008 PHYSICAL REVIEW LETTERS In this Letter we introduce a new stochastic Lagrangian model for flow and transport in porous media that allows one to separate advective mixing from its diffusive counterpart The model posits that fluid flow in homogeneous porous media is governed by a combination of the continuity equation 1 and a stochastic Langevin flow equation which is obtained by adding white noise fluctuations to the macroscopic flow equation 2 In this mesoscale formulation the noise represents the subgrid variability of the flow velocity and the combined effects of the simplifying assumptions leading to 2 and accounts for deviations from the smooth flow paths predicted by the Darcy scale continuum flow equations 1 and 2 This is conceptually analogous to the role played by the noise in Brownian motion models of diffusion and in fluctuating hydrodynamics 8 The resulting stochastic flow and transport equations can be solved by a variety of methods Here we solve them with smoothed particles hydrodynamics SPH a Lagrangian numerical algorithm that has been successfully applied to both deterministic 9 10 and stochastic 11 transport problems In SPH the fluid is represented with M particles and the stochastic Langevin flow equation is q rPi dXi dUi Ui g i Ui jhUi ij i dt dt 4 where the subscript i i 1 M indicates the random position X velocity U and pressure P of the ith particle Reynolds decomposition is used to represent a random of quantity A X U P as the sum A hAi A its ensemble mean hAi and random fluctuations about the The components of the white noise mean A 1 d T satisfy h l t i 0 l 1 d where d is the system s dimensionality and 2 t t0 l m 0 5 h l t m t i ll 0 otherwise For simplicity and without any loss of generality we set the constants ll for all l This is consistent with taking the dispersivity to be a scalar The variance of the ith particle s position can be related to the dispersivity via Einstein s relationship dhX l i t X l i t i 2 jhUi i1 j t 1 dt lim 6 where hUi i1 is the steady average velocity of the ith particle In the stochastic simulations presented here the