Foams are used in reservoir engineering for enhanced oil recovery, CO2 sequestration and environmental remediation of aquifers and soils. One of the main mechanisms for foam generation at steady state is Roof snap-off. In some cases, mechanistic models of Roof snap off, based on observations from 2D micromodels are used for reservoir simulation. The main problem with these experiments is in their 2D nature. Two-phase flow within a 2D medium requires that the fluids paths cross and compete for pore occupancy. This virtually ensures fluctuating pore occupancy and therefore puts into question the applicability of 2D mechanistic models for steady state foam generation in 3D media. Two-phase flow in a micromodel is analyzed with a lattice percolation model in order to determine under what conditions steady two-phase flow can be achieved. The gas network is established with bond percolation and liquid is allowed to flow across the sample with the help of liquid bridges. These liquid bridges enable the liquid to cross gas-occupied pore throats without snap-off. The calculated attribute for the gas and liquid networks is equivalent resistance, ΔP/Q. For this a new unit for hydraulic resistivity was used and is equal to the fluid viscosity divided by the pore radius to the third power, H = μ/R3. The equivalent resistance of the gas and liquid networks are found by applying rules from linear circuits of electrical resistances. Solutions for the equivalent resistivity of the gas network are calculated with the node elimination method and Kirchhoff’s solution for a random network of resistances. The liquid network’s conductivity is calculated as the sum of path resistances in parallel and is a theoretical maximum. The gas and liquid conductivity of nine pre-existing 16x16 networks from Holstvoogd(2020) are reevaluated and, in addition, twelve new samples of size 32x32 are evaluated. Functionally, the model’s behavior is as follows: gas conductivity is inversely proportional and liquid conductivity is proportional to the occupation threshold. It is found that the gas conductivity is a function of tortuosity and number of parallel flow loops. Conductivity decreases with increased tortuosity and increases with number of parallel flow paths. The ratio of liquid and gas conductivity for the twelve 32x32 models is calculated. When adjusted for gas viscosities of supercritical CO2 and Nitrogen gas it is found that it is in the order of 10-3 to 10-4. Therefore, it has been determined that it is practically impossible to achieve steady two-phase flow without fluctuating pore occupancy.

Compositional simulation is computationally intensive for high-fidelity models due to thermodynamic equilibrium relations and the coupling of flow, transport and mass transfer. In this report, two methods for accelerated compositional simulation are outlined and demonstrated for a gas vaporization problem. The first method uses a proxy model that reduces the number of components and the second method reduces the number of grid blocks (i.e. upscaling). Both methods are implemented within the operator-based linearization framework of the Delft Advanced Research Terra Simulator. Lebesgue integration is applied in the loss function of a neural network allowing the neural network to discover the operator space of the reference model in reduced dimensions. Training is carried out for a one-dimensional homogeneous reservoir and minimizes the misfit of the leading and trailing shocks of a compressible pseudo-binary model with respect to observations of the reference model. The operator space of the pseudo-binary model is initially approximated with the method of multiscale reconstruction of physics, a numerical representation of the method of characteristics. Training is carried out in a two-stage transfer learning scheme to increase computational efficiency. In the first stage, neural networks are trained to approximate the analytical reconstruction. In the second stage, a solver is embedded in the loss function of the neural network and the forward solution is used to calculate the Lebesgue integral. The transfer training scheme minimizes the misfit of the leading and trailing shocks for 10 discrete time steps in a one-dimensional homogeneous reservoir. The misfit of the trained model shows a significant improvement in the location of the trailing shock and a modest improvement in the estimation of the leading shock. The trained proxy is applied to the top and bottom 15 layers of the SPE10 model and the estimation of the first and last breakthrough is assessed in conjunction with the error of the phase-state classification. The phase-state classification is significantly improved through time which is also expressed in improvements of the estimation of breakthrough times. The average difference in breakthrough time for the trained and untrained models with respect to the reference model is 293days versus 570days for the trailing shocks and 15days versus 16days for the leading shock. The established training framework enables the development of proxies with increased complexity. Rigorous upscaling defines the upscaled operator space with dynamic and non-equilibrium thermodynamic upscaling functions. These functions combined, define the upscaled operator space for the three-dimensional compositional space and are inferred from data points gathered from a limited, characteristic portion of the full-size model. Gathered data points are interpreted with an interpolation function or neural networks to construct structured OBL meshes for implementation within DARTS. This upscaled operator space can effectively be used for different boundary conditions without reevaluating the upscaling functions.