This dissertation focuses on the application of Monte Carlo (MC) molecular simulation and Molecular Dynamics (MD) in modeling thermodynamics and flow of subsurface reservoir fluids. At first, MC molecular simulation is proposed as a promising method to replace correlations and equations of state in subsurface flow simulators. In order to accelerate MC simulations, a set of early rejection schemes (conservative, hybrid, and non-conservative) in addition to extrapolation methods through reweighting and reconstruction of pre-generated MC Markov chains were developed. Furthermore, an extensive study was conducted to investigate sorption and transport processes of methane, carbon dioxide, water, and their mixtures in the inorganic part of shale using both MC and MD simulations. These simulations covered a wide range of thermodynamic conditions, pore sizes, and fluid compositions shedding light on several interesting findings. For example, the possibility to have more carbon dioxide adsorbed with more preadsorbed water concentrations at relatively large basal spaces. The dissertation is divided into four chapters. The first chapter corresponds to the introductory part where a brief background about molecular simulation and motivations are given. The second chapter is devoted to discuss the theoretical aspects and methodology of the proposed MC speeding up techniques in addition to the corresponding results leading to the successful multi-scale simulation of the compressible single-phase flow scenario. In chapter 3, the results regarding our extensive study on shale gas at laboratory conditions are reported. At the fourth and last chapter, we end the dissertation with few concluding remarks highlighting the key findings and summarizing the future directions.

Third phase formation is an undesirable phenomenon during the PUREX process, which is a continuous liquid-liquid extraction approach for the reprocessing of uranium and plutonium from spent nuclear fuel. When third phase formation occurs, the organic extraction solution splits into two layers. The light upper layer, which is commonly named the light organic phase, contains a lower concentration of metal ions, tri-n-butyl phosphate (TBP) and nitric acids but is rich in the organic diluent. The heavy lower layer, which is commonly named the third phase, contains high concentrations of metal ions, TBP and nitric acids. As the third phase contains high concentrations of the uranium and plutonium complexes it can thus cause processing and safety concerns. Therefore, a comprehensive understanding of the mechanism of third phase formation is needed so as to improve the PUREX flowsheet. To investigate third phase formation through molecular simulations, one should first obtain reliable molecular models. A refined model for TBP, which uses a new set of partial charges generated from our density functional theory calculations, was proposed in this study. To compare its performance with other available TBP models, molecular dynamics simulations were conducted to calculate the thermodynamic properties, transport properties and the microscopic structures of liquid TBP, TBP/water mixtures and TBP/n-alkane mixtures. To our knowledge, it is only TBP model that has been validated to show a good prediction of the microscopic structure of systems that consist of both hydrophobic and hydrophilic species. This thesis also presents evidence that the light-organic/third phase transition in the TBP/n-dodecane/HNO3/H2O systems, which is relevant to the PUREX process, is an unusual transition between two isotropic, bi-continuous micro-emulsion phases. The light-organic /third phase coexistence was first observed using Gibbs Ensemble Monte Carlo (GEMC) simulations and then validated through Gibbs free energy calculations. Snapshots from the simulations as well as the cluster analysis of the light organic and third phases reveal structures akin to bi-continuous micro-emulsion phases, where the polar species reside within a mesh whose surface consists of amphiphilic TBP molecules. The non-polar n-dodecane molecules are outside this mesh. The large-scale structural differences between the two phases lie solely in the dimensions of the mesh. To our knowledge, the observation of the light-organic/third phase coexistence through simulation approaches and a phase transition of this nature have not previously been reported. Finally, this thesis presents evidence that the microscopic structure of the light organic phase of the Zr(IV)/TBP/n-octane/HNO3/H2O system, which is also related to the PUREX process, is different from that of the common hypothesis, where such system is consisted of large ellipsoidal reverse micelles. Snapshots from simulations, hydrogen bonding analysis and cluster analysis showed that the Zr4+, nitrate, TBP and H2O form extended aggregated networks. Thus, as above, we observe a bi-continuous structure but this time with embedded local clusters centred around the Zr4+ ions. The local clusters were found to consist primarily of Zr(NO3)4·3TBP complexes. This finding provides a new view of the structure of the Zr(IV)/TBP/n-octane/HNO3/H2O system.

Homogeneous crystal nucleation in binary hard sphere mixtures is an active area of research for last two decades. Although Classical nucleation theory (CNT) gives a qualitative picture, it fails at high super saturations because of the following reasons. CNT assumes that the cluster formed is spherical in shape, its properties can be modeled using the bulk properties of the stable solid phase and the interfacial free energy γ between the nucleus and the surrounding metastable fluid is equal to the planar surface tension between two phases at coexistence. These assumptions get increasingly tenuous at higher degrees of super saturations where the critical nucleus formed is microscopic in size leading to breakdown in the predictions of CNT. In addition direct experimental observation of critical nucleus is very difficult because, 1. Critical nucleus is microscopic in size, consisting of few hundreds of particles. 2. Formation of critical cluster is very rare (typically of the order of 101– 106nuclei/cm3/s) 3. Its life time is very short (it either rapidly grows to form a solid phase or melts back to fluid) In these circumstances molecular simulations are an attractive tool to study the crystal nucleation, because in these simulations microscopic size critical nucleus properties can be calculated. However, brute force molecular dynamic (MD) simulation techniques to study the homogeneous crystal nucleation is currently not feasible due to long times involved. Hence, an indirect approach is needed. In this work, Monte Carlo Abstract v (MC) molecular simulation techniques are used to calculate free energy barrier height during the crystal nucleation. Phase behavior of Binary hard sphere mixtures with varying ratios of smaller diameter to larger diameter (α) is very similar to that of binary organic liquids. By studying the crystal nucleation in hard sphere system, the physics behind the nucleation for binary organic liquids can be understood. This is the key motivation to study the homogeneous crystal nucleation in binary hard sphere mixtures using MC simulations. Simulations were done using umbrella sampling in combination with local bond order analysis for the identification of crystal nuclei and to compute shape and height of nucleation barrier. Parallel tempering scheme of Geyer and Thomson was utilized to sample phase space more efficiently. Parallel tempering technique was implemented using Message passing interface (MPI) libraries. By using all the above Monte-Carlo simulation techniques, nucleation barrier was calculated during crystallization of binary hard sphere mixtures under the moderate degrees of super cooling in Isothermal-Isobaric semi grand ensembles. Crystal nucleation in binary hard sphere mixtures has been studied for size ratios α = 0.85, 0.42 and 0.43. For α=0.85, phase diagram contains eutectic point. In this system, the effect of eutectic composition on the nucleation barrier height was observed, by calculating nucleation barriers at various fluid mixture compositions keeping Laplace pressure constant. It is observed that as the fluid mixture composition move towards the eutectic point, free energy barrier height, surface tension and critical cluster sizes are increased and the nucleation rate is drastically decreased by a factor of 10-31. Thus the difficulty of homogenous crystal nucleation increases near the eutectic point. For α=0.42 and 0.43 in the hard sphere system, compound solids such as AB and AB2 solids are stable respectively. In these systems crystal nucleation study was done to observe the compound solid formation. It is observed that in these systems crystallization kinetics are very slow and more advanced simulation techniques need to be developed in order to study crystal nucleation.

Semicrystalline polymers are an important class of materials for their prevalence in today's markets and their desirable properties. These properties depend on the early stages of the polymer crystallization process where a crystal nucleates from the polymer melt. This nucleation process is conventionally understood via an extension of Classical Nucleation Theory to polymers (CNTP). However, recent experimental and simulation evidence points to nucleation mechanisms that do not agree with the predictions of CNTP. Specifically, these experiments suggest a previously unrecognized role of nematic phases in mediating the melt"“crystal transtion. To explain these observations, several new theories of nucleation alternate to CNTP have emerged in the literature, all of which suggest specific modifications to the free energy landscape (FEL) near-equilibrium. To address these theoretical controversies, this dissertation aimed to study the equilibrium phase behavior of polymers via Monte Carlo (MC) simulations. Simulating equilibrium phase behavior of polymer melts is not a trivial task due to the large free energy barriers involved. Throughout this research, we employed a combination of strategies to speed up these molecular simulations. First, we employed a domain decomposition to divide the simulation box into multiple independent simulations that execute independent MC trajectories in parallel. The novel GPU-accelerated MC algorithm successfully and accurately simulated the phase behavior of bead spring chains. Additionally, it sped up MC simulations of Lennard Jones chains by up to 10 times. In its current form, the GPU-accelerated algorithm did not achieve significant speedups to improve outcomes of simulating large polymer melts with detailed potentials. We recommended various strategies to improving the current algorithm. This reality motivated the use of biased MC simulations to study the phase behavior of polymers more expediently without the need for GPU acceleration. Specifically, the latter part of the Dissertation employed Wang Landau MC (WLMC) simulations to build phase diagrams and expanded ensemble density of states (EXEDOS) simulations to construct FELs. Phase diagrams from WLMC simulations divided volume-temperature space into melt, nematic and crystal phases. Then, FELs from EXEDOS simulations at equilibrium provided direct access to the relative stability and minimum free energy paths between coexistant states. By employing a two-dimensional EXEDOS sampling in both crystal and nematic order for hard bead semiflexible oligomers with a stepwise bending stiffness, we built FELs that show that the crystalline transition cooperatively and simultaneously formed crystal and nematic order. This nucleation mechanism was not in agreement with predictions from CNTP or newer theoretical formulations. To investigate the sensitivity of the phase behavior to the employed polymer model, we then employed WLMC simulations to build phase diagrams for a number of different polymer models to ascertain their impact on the resulting nucleation mechanism. We found that the phase behavior was sensitive to the form of the bending stiffness potential used. Chains with a stepwise bending stiffness yielded the previously mentioned cooperative and simultaneous crystal and nematic ordering. In contrast, chains with a harmonic bending stiffness potential crystallized via a two-step nucleation process, first forming a nematic phase that nucleates the crystal. The latter nucleation mechanism was in line with predictions from new theories of nucleation that incorporate the nematic phase as a precursor. Furthermore, we found that it is important to correct for excluded volume differences when comparing chains with soft and hard beads or chains with differing bending stiffnesses.