Pore-network models have become a critical tool in the study of fluid flow in geo-energy researches over the last few decades, and the accuracy of pore-network modeling results highly depends on the extraction of pore networks. Traditional methods of pore-network extraction are based on pixels and require images with high quality. Here, a pixel-free method called the flashlight search medial axis algorithm is proposed for pore-network extraction in a continuous space. The search domain in a two-dimensional space is a line, whereas a surface domain is searched in a three-dimensional scenario. Thus, the algorithm follows the dimensionality reduction idea; the medial axis can be identified using only a few points instead of calculating every point in the void space. In this way, computational complexity of this method is greatly reduced compared to that of traditional pixel-based extraction methods, thus enabling large-scale pore-network extraction. Based on cases featuring two- and three-dimensional porous media, the algorithm performs well regardless of the topological structure of the pore network or the positions of the pore and throat centers. This algorithm can also be used to examine both closed-boundary and open-boundary cases. Finally, this algorithm can identify the medial axis accurately, which is of great significance in the study of geo-energy.

Spontaneous imbibition is a crucial process for oil recovery from fractured and unconventional reservoirs. Herein, with the assumption of capillaries being independent, a new mathematical model for spontaneous imbibition is proposed and solved using a numerical method. The simulated results show that the wetting phase preferentially enters smaller capillaries where the advancement velocity is higher than that in larger ones, while the non-wetting phase can be displaced out in the larger capillaries. In addition, the effect of fluid viscosity ratio on counter-current imbibition is analyzed. The results show that imbibition velocity becomes higher with the increase in the viscosity ratio. When the viscosity of the non-wetting phase is larger than that of the wetting phase, the end pressure gradually increases as the imbibition front advances. In contrast, when the viscosity of the non-wetting phase is less than that of the wetting phase, the end pressure decreases with the infiltration. With a higher viscosity ratio of non-wetting and wetting phase, the heterogeneity of the interface advancement among different capillaries increases.

In the last decades, shale gas development has relieved the global energy crisis and slowed global warming problems. The water bridge plays an important role in the process of shale gas diffusion, but the stability of the water bridge in the shale nanochannel has not been revealed. In this work, the molecular dynamics method is applied to study the interaction between shale gas and water bridge, and the stability can be tested accordingly. CO₂ can diffuse into the liquid H₂O phase, but CH₄ only diffuses at the boundary of the H₂O phase. Due to the polarity of H₂O molecules, the water bridge presents the wetting condition according to model snapshots and one-dimensional analyses, but the main body of the water bridge in the two-dimensional contour shows the non-wetting condition, which is reasonable. Due to the effect of the molecular polarity, CO₂ prefers to diffuse into kerogen matrixes and the bulk phase of water bridge. In the bulk of the water bridge, where the interaction is weaker, CO₂ has a lower energy state, implies that it has a good solubility in the liquid H₂O phase. Higher temperature does not facilitate the diffusion of CO₂ molecules, and higher pressure brings more CO₂ molecules and enhances the solubility of CO₂ in the H₂O phase, in addition, a larger ratio of CO₂ increases its content, which does the same effects with higher pressures. The stability of the water bridge is disturbed by diffused CO₂, and its waist is the weakest position by the potential energy distribution.

In order to make full use of the information provided in the physical reservoirs, including the production history and environmental conditions, the whole life cycle of reservoir discovery and recovery should be considered when mapping in the virtual space. A new concept of reservoir digital twin and the exploratory multi-scale framework is proposed in this paper, covering a wide range of engineering processes related with the reservoirs, including the drainage, sorption and phase change in the reservoirs, as well as extended processes like injection, transportation and on-field processing. The mathematical tool package for constructing the numerical description in the digital space for various engineering processes in the physical space is equipped with certain advanced models and algorithms developed by ourselves. For a macroscopic flow problem, we can model it either in the Navier-Stokes scheme, suitable for the injection, transportation and oil-water separation processes, or in the Darcy scheme, suitable for the drainage and sorption processes. Lattice Boltzmann method can also be developed as a special discretization of the Navier-Stokes scheme, which is easy to be coupled with multiple distributions, for example, temperature field, and a rigorous Chapman-Enskog expansion is performed to show the equivalence between the lattice Bhatnagar-Gross-Krook formulation and the corresponding Navier-Stokes equations and other macroscopic models. Based on the mathematical toolpackage, for various practical applications in petroleum engineering related with reservoirs, we can always find the suitable numerical tools to construct a digital twin to simulate the operations, design the facilities and optimize the processes.

Compositional multiphase flow in subsurface porous media is becoming increasingly attractive due to issues related with enhanced oil recovery, CO ${}_{\mathrm{2}}$ sequestration and the urgent need for development in unconventional oil/gas reservoirs. One key effort to construct the mathematical model governing the compositional flow is to determine the phase compositions of the fluid mixture, and then calculate other related physical properties. In this paper, recent progress on phase equilibrium calculations in unconventional reservoirs has been reviewed and concluded with authors’ own analysis, especially focusing on the special mechanisms involved. Phase equilibrium calculation is the main approach to investigate phase behaviors, which could be conducted using different variable specifications, such as the NPT flash and NVT flash. Recently, diffuse interface models, which have been proved to possess a high consistency with thermodynamic laws, have been introduced in the phase equilibrium calculation, incorporating the realistic equation of state (EOS), e.g. Peng-Robinson EOS. In the NVT flash, the Helmholtz free energy is minimized instead of the Gibbs free energy used in NPT flash, and this thermodynamic state function is decomposed into two terms using the convex-concave splitting technique. A semi-implicit numerical scheme is applied to the dynamic model, which ensures the thermodynamic stability and then preserves the fast convergence property. A positive definite coefficient matrix is designed to meet the Onsager reciprocal principle so as to keep the entropy increasing property in the presence of capillary pressure, which is required by the second law of thermodynamics. The robustness of the proposed algorithm is demonstrated by using two numerical examples, one of which has up to seven components. In the complex fluid mixture, special phenomena could be captured from the global minimum of tangent plane distance functions and the phase envelope. It can be found that the boundary between the single-phase and vapor-liquid two phase regions shifts in the presence of capillary pressure, and then the area of each region changes accordingly. Furthermore, the effect of the nanopore size distribution on the phase behavior has been analyzed and a multi-scale scheme is presented based on literature reviews. Fluid properties including swelling factor, criticality, bubble point and volumetrics have been investigated thoroughly by comparing with the bulk fluid flow in a free channel.

In this paper we review the Dynamic Van der Waals theory, which is a recent developed method to study phase separation and transition process in multiphase flow. Gradient contributions are included in the entropy and energy functions, and it’s particularly useful and non-trivial if we consider problems with temperature change. Using this theory, we can simulate that, a droplet in an equilibrium liquid will be attracted to the heated wall(s) which was initially wetted, which is the main cause of the famous hydrodynamic phenomena-Leidonfrost Phenomena. After more than ten years development, this theory has been widely used to study the fluid flow in vaporing and boiling process, e.g., droplet motion. Furthermore, this theory has been combined with phase field model, which could be extended to solid-liquid phase transition. There has also been researches about constructing LBM scheme to extend to the Dynamic Van der Waals theory, using Chapman-Enskog analyze. In all, due to its rigorous thermodynamic derivation, this theory has now become the fundamental theoretical basis in the heated multiphase flow.