Viscous fingering occurs when a less viscous fluid is injected into a rock matrix saturated with a more viscous fluid. Our past research using the Rothman-Keller (RK) color gradient Lattice Boltzmann Method (LBM) for immiscible two phase flow has allowed us to study viscous fingering morphology and the complex saturation phase space as a function of the fluid’s properties (wettability of the injected fluid and viscosity ratio). In this past work, we found that the primary factor affecting the saturation at breakthrough – when the injected fluid has passed through the entire model – was the viscosity ratio, and the secondary effect was the wettability. Here, we present an extension of our LBM model to enable convection-diffusion to be simulated, thereby allowing us to vary the viscosity of the injected fluid, and mimicking the practice in Enhanced Oil Recovery (EOR) using polymer additives after breakthrough as a means of increasing the viscosity ratio and thus the eventual oil yield. The basic RK multiphase LBM models two fluid number densities moving and colliding on a discrete lattice, where a second collision term is used to model cohesion within each fluid, and contains an extra “recoloring step” to ensure fluid segregation. Here, we model an additional number density representing the concentration of a polymer additive, which affects the viscosity of the injected fluid. The Peclet number – rate of advection to diffusion of the polymer solution – is used to set the diffusion coefficient of the polymer concentration number density and hence, the relaxation time in the LBM for the polymer diffusion process. We present tests to demonstrate the method in which we increase the polymer concentration of the injected fluid after a given time and study the effect on the viscous fingering morphology and saturation evolution. This work demonstrates that the RK color gradient multiphase LBM can be used to study complex viscous fingering behavior associated with injection of water with polymer additives, which can have major scientific and practical significance.

Recent advances in modeling Rayleigh-Benard convection have demonstrated the existence of turbulent superstructures, whose life and morphology largely varies with Rayleigh (Ra) and Prandtl (Pr) numbers. These structures appear as a two scale phenomena, where small scale rolls organize in larger convection cells, and can be modelled only in 3D on a simulation box characterized by a very large (>10) width/height (W/L) ratio, and sufficiently refined to resolve the boundary layer up to Ra = 108 (>100 divisions in height) and to Ra = 1010 (>200 divisions). To achieve this goal, we use our own 3D Parallel Python implementation of the Lattice Boltzmann Method, tested to run with linear efficiency on thousands of cores. We show the dependency of horizontal fluctuations of RMS of temperature and vertical velocity in the middle of the box and illustrate how the superstructures emerge for W/L ratios of Terrestrial Planets and Super Earths, and quantify the duration of these superstructures and the likely implications for the evolution of their surface features. The effect of the P/T dependent viscosity and thermal conductivity is finally discussed.

Lithospheric yield stress is a key parameter in controlling tectonic processes. We calculate yield stress for a range of conditions appropriate to the Archean Earth, including hotter mantle potential temperatures and a range of Moho temperatures using 2D high resolution numerical geodynamic modelling techniques. This range of conditions are evaluated for generating felsic, tonalite-trondhjemite-granodiorite (TTG), crust with the results bench marked against the preserved rock record. The model results indicate that lithospheric yield stress slightly lower than the present-day Earth values (i.e. < 100 MPa) generates TTG melt volumes similar to those preserved in the rock record. In particular, large volumes of TTG melts form in the tails of lithospheric drips. Melting occurs profusely within the thinner portions of the drips as these regions are more ef�ciently heated by the enclosing hotter mantle. In contrast, only limited melting occurs in regions of thickened crust, in part because the weaker lithosphere cannot sustain crustal thickening for long time periods, resulting in its removal through drips. Our models highlight the dominance of non-plate tectonic mechanisms in producing TTGs under the conditions that operated on the hotter Archean Earth.

Cratons are stable parts of the Earth’s continental lithosphere that have remained largely undeformed for several billion years. These consist of crustal granite-greenstone terrains coupled to roots of strong, buoyant cratonic lithospheric mantle that extend up to several hundreds of kms depth. Due to their stability, cratons preserve a record of the tectonics and the thermal evolution of the mantle in the early Earth. These observations suggest that the highly viscous (strong) character of cratons hampered the viability of early Earth tectonics, thereby affecting mantle convection patterns and cooling. In this study, we investigate the controls of stiff cratons on the initiation of subduction and mantle thermal evolution on the early Earth. Using numerical models, we simulate the effects of strong and buoyant cratons on mantle convection. We vary a set of parameters including (i) width and thickness of cratons, and (ii) viscosity ratio between cratonic lithosphere and cratonic crust. We test initial conditions varying the number of cratons, which is unconstrained for early Earth and associated it to mantle cooling rates. Our preliminary results show that the mantle cooling rate decreases with increasing number of cratons. Because mantle cooling rates affect the early Earth transition from a basaltic drip regime to initiation of subduction, we show that the craton coverage on the early Earth controls the time of onset of plate tectonics. Furthermore, we observe that cratons will remain separate or combine depending on the convective cell size, which is function of mantle cooling.