A realistic deep low-enthalpy geothermal reservoir based on real data with high detail and complicated sedimentary structure is utilized to perform sensitivity analyses of the geological features influencing reservoir properties. We perform simulations using the Delft Advanced Research Terra Simulator (DARTS). Compelling numerical performance of DARTS makes it suitable for handling a large ensemble of models including efficient sensitivity and uncertainty analyses. The major finding is that shale facies, generally ignored in hydrocarbon reservoir simulations, can significantly extend the predictive lifetime of geothermal reservoirs exploited by deep well doublets. It is important to accurately account for the shale facies in the simulation, though with an additional computational overhead. The overburden layers can improve doublet performance, but the impact depends on reservoir heterogeneity. In addition, heterogeneity will also divert the flow path with even a minor shift in the well placement. The discharge rate, an essential parameter of geothermal operation strategy, inversely corresponds to the doublet lifetime but positively correlates with the energy production for studied parameter ranges. Low sensitivity of doublet lifetime to vertical-horizontal permeability ratio and permeability-porosity correlation is observed. All these systematic findings for a realistic geothermal field with characterization at unprecedented level of detail can help to provide a general guideline for forward simulation and farther improve the profitability of geothermal energy production in realistic deep geothermal reservoirs through computer-assisted modeling and optimization.

@en

A physics-based data-driven model is proposed for forecasting of subsurface energy production. The model fully relies on production data and does not require any in-depth knowledge of reservoir geology or governing physics. In the proposed approach, we use the Delft Advanced Reservoir Terra Simulator (DARTS) as a workhorse for data-driven simulation. DARTS uses an operator-based linearization technique that exploits an abstract interpretation of physics benefiting computational performance. The physics-based data-driven model is trained to fit data increasing the fidelity of the model forecast and reflecting significant changes in reservoir dynamics or physics over its history. The model is examined and validated for both synthetic and real field production data. We demonstrate that the developed approach is capable of providing accurate and reliable production forecast on a daily basis, even if the exact geological information is not available.

@en

Accurate prediction of temperature and pressure distribution is essential for geothermal reservoir exploitation with cold water re-injection. Depending on our knowledge about the heterogeneous structure of the subsurface, the reservoir development scheme can be optimized and the overall lifetime of the geothermal field can be extended. In this study, we present Delft Advanced Research Terra Simulator (DARTS), which provides fast and accurate energy production evaluation for geothermal applications. This simulation framework is suitable for uncertainty analysis with a large ensemble of models. In DARTS, we select the molar formulation with pressure and enthalpy as primary variables. Besides, the fully-coupled fully-implicit two-point flux approximation on unstructured grids is implemented to solve the mass and energy conservation equations. For the nonlinear solution, we employ the recently developed Operator-Based Linearization (OBL) approach. In our work, DARTS is compared with the state-of-the-art simulation frameworks using a set of benchmark tests. We demonstrate that DARTS achieves a good match for both low- and high-enthalpy conditions in comparison to other simulators. At the same time, DARTS provides high performance and flexibility of the code due to the OBL approach, which makes it particularly useful for uncertainty quantification in processes involving complex physics.

@en

Numerical modelling of multiphase multicomponent flow coupled with mass and energy transport in porous media is crucially important for many applications including oil recovery, carbon storage and geothermal. To deliver robust simulation results, a fully or adaptive implicit method is usually employed, creating a highly nonlinear system of equations. It is then solved with the Newton-Raphson method, which requires a linearization procedure to assemble a Jacobian matrix. Operator Based Linearization (OBL) approach allows detaching property computations from the linearization stage by using piece-wise multilinear approximations of state-dependent operators related to complex physics. The values of operators used for interpolation are computed adaptively in the parameter-space domain, which is uniformly discretized with the desired accuracy. As the result, the simulation performance does not depend on the cost of property computations, making it possible to use expensive equation-of-state formulations (e.g., fugacity-activity thermodynamic models) or even black-box chemical packages (e.g., PHREEQC) for an accurate representation of governing physics without penalizing runtime. On the other hand, the implementation of the simulation framework is significantly simplified, which allows improving the simulation performance further by executing the complete simulation loop on GPU architecture. The integrated open-source framework Delft Advanced Research Terra Simulator (DARTS) is built around the OBL concept and provides a flexible, modular and computationally efficient modelling package. In this work, we evaluate the computational performance of DARTS for various subsurface applications of practical interests on both CPU and GPU platforms. We provide a detailed performance comparison of particular workflow pieces composing Jacobian assembly and linear system solution, including both stages of Constrained Pressure Residual solver.

@en

The efficient development of a geothermal field can be largely affected by the inherent geological and physical uncertainties. Besides, the uncertain operational and economic parameters can also impact the profit of a project. Systematic uncertainty quantification involving these parameters helps to determine the probability of concerning outputs. In this study, a low-enthalpy geothermal reservoir with strong heterogeneity, located in the West Netherlands Basin, is selected as the research area. Detailed geological model is constructed based on various static data including seismic and log interpretation. However, significant uncertainties still exist in definition of the model parameters, mainly reservoir permeability and porosity. Besides, the fluid properties have not been sampled in this field and can vary in the range between brackish to highly saline water. Also, the heat price and operational investment fluctuate with time and add up to uncertainty. Taking all interested parameters into consideration, the Monte Carlo method is utilized to select specific input data set. The forward simulations are powered by the GPU version of Delft Advance Research Terra Simulator (DARTS), which provides efficient simulation capabilities for geothermal applications. Through this investigation, a wide range of production temperature has been observed due to the uncertainty of the input parameters.@en

The continuous progress of reservoir monitoring technology provides encouraging capacities to reduceuncertainties in the subsurface characterization and to mitigate risks in field development applying thereservoir simulation approach. However, it is always challenging to take full advantage of the observationdata, since an accurate representation of strong heterogeneities requires a high-resolution grid. Most ofthe discretization methods cannot handle full tensor permeability, and high nonlinearity introduced bycomplex physical process drastically reduces simulation efficiency. In this work, we develop an advancedparallel framework for reservoir simulation with the implementation of state of the art discretizationand linearization methods. We apply the multipoint flux approximation (MPFA) method to handle thefull tensor permeability in unstructured grids. To keep the fidelity of the geological model and improvecomputational efficiency, we use massively parallel computations via Message Passing Interface (MPI).Complex subsurface physics is described by mass-based formulations making the framework flexible forgeneral-purpose reservoir simulation. However, the representation of phase behavior introduces additionalworkload when compared with the phase-based formulations in the traditional approach. Here, we apply theOperator-Based Linearization (OBL) approach which not only overcomes this drawback but also turns it toan advantage. In this method, the conservation equations are described in an operator form. By constructinga library of tabulated operators, the repeated work spent on complex phase behavior and property evaluationcan be significantly reduced. We benchmark the parallel framework with analytical solutions under single-phase flow and multiphase flow. The results demonstrate that the parallel framework provides accuratesimulation results for structured and unstructured grids. We validate that MPFA implemented in our parallelframework converges to real solutions when the permeability is a full tensor. Besides, several realisticcases have been rigorously tested confirming high computational capacity, efficiency, and accuracy of theadvanced massively parallel framework for general-purpose reservoir simulation. With the implementationof MPFA and OBL approaches, the parallel framework is fully equipped for the simulation of problemswith full tensor permeability, high-heterogeneities, and complex physical processes.

@en

A physics-based data-driven model is proposed in this study for the forecasting of secondary oil recovery. The model fully relies on production data and does not directly requires any in-depth knowledge of the reservoir geology or governing physics. In the proposed approach, we utilise Delft Advanced Reservoir Terra Simulator (DARTS) as a base for data-driven simulation. DARTS uses an Operator-Based Linearization technique which exploits an abstract interpretation of physics benefiting computational performance for a forward simulation. The proposed strategy was evaluated first on the two synthetic data ensembles and showed good prediction accuracy for a significantly reduced model size. Besides, the data-driven proxy methodology was compared with an advanced flow-based upscaling technique and demonstrated an improved accuracy for both ensembles. Besides, the proposed data-driven approach was examined on two realistic data sets. For the first case, the methodology demonstrates advanced predictive performance for training based on synthetic data generated from a high-fidelity simulation model with imposed random noise. To check the robustness of the proposed methodology, the control parameters for a forecast period were significantly changed in comparison to the training period. The data-driven model still manages to predict the forecast production quite close to the reference high-fidelity results. However, the training performed on another data set based on historical production from a real brownfield was not fully successful. We relate a bigger error in both training and forecast period for this model to poor data quality. The training procedure for this model led to a moderate accuracy in history matching for a long production period, where general production trends have resembled true data and water breakthrough time was restored in nearly all wells. However, there are still periods of poor accuracy, especially where shark peaks and falls are experienced.

@en

Numerical simulation is based on space and time discretization providing a trade-off between accuracy and computational performance. The operator-based Linearization (OBL) approach introduces an additional discretization domain for the physical description of fluid and rock. Delft Advanced Research Terra Simulator (DARTS) is a general-purpose reservoir simulation platform entirely built around the OBL approach. DARTS provides unique flexibility capabilities, allowing to customize physical description and control simulation process using the high-level Python programming language. Fully Implicit Method provides unconditional simulation stability, while the partial derivatives required for Jacobian assembly are computed analytically through OBL. At the same time, DARTS ensures exceptional simulation performance addressing it at three levels. Inexpensive linearization combined with a reduction in the nonlinearity is controlled by OBL discretization resolution form algorithmic level. Efficient C++ backend for critical kernels and state-of-the-art linear solver with two-stage CPR preconditioner composes the software level. Finally, the ability to perform the entire simulation loop on the GPU platform owing to the OBL approach constitutes the hardware level. The DARTS framework has already served as a platform for several academic and industrial research projects for different geo-energy applications including geothermal, CO2 sequestration and petroleum.@en

Various novel computing architectures, like massively parallel and multi-core, as well as computing accelerators, like GPUs or TPUs, keep regularly expanding. In order to exploit the benefits of these architectures to the full extent and speed up reservoir simulation, the source code has to be inevitably rewritten, sometimes almost completely. We demonstrate how to extract complex physics-related computations from the main simulation loop, leaving only an algebraic multilinear interpolation kernel instead. In combination with linear solvers, which usually have made available soon once the new architecture is introduced, the approach accommodates execution of the entire nonlinear loop on the latest hardware and computational architectures. We describe the integrated simulation framework built on top of this technique and show the applicability of the approach to various challenging physical and chemical problems. All simulation engines along with linear solvers, well controls, interpolation engines, and state operator evaluators are implemented in C++11 and exposed into Python coupling the flexibility of the script language with the performance of C++.

@en

Petrographic description is one of the primary methods by which geoscientists develop an understanding of geologic systems and constrain source, reservoir, and seal characteristics. Geologic facies or rock types are routinely interpreted from petrographic features in thin sections, but these methods are time consuming and tedious. The interpretations can be heavily biased. Variations in experience base, microscope quality, and scale of observation can lead to inconsistent interpretations, absorbing costly time to achieve alignment. With advances in computing, millions of images can be analyzed in seconds through signal-processing technologies. With this capability, machine-learning methods can support geologic pattern identification to classify geologic features. Combined, automated image analysis and classification though machine learning can significantly reduce analysis time, interpreter bias and inconsistencies. These methods make it possible to share "expert knowledge" with nonexpert users. Image processing, segmentation and classification are key workflows for translating geologic features into discrete representations that can be used for computational modeling. Here, we present two examples where advanced image analysis and machine learning are used to predict geologic and petrophysical properties from optical microscopy thin-section images. The first example uses multiple machine-learning algorithms and pore geometries as a means for predicting and classifying rock properties, such as lithofacies, reservoir zone, porosity and permeability. The second example focuses on training and using convolutional neural networks to classify and predict Dunham textures (mudstone, wackestone, packstone and grainstone) from carbonate thin sections. Data for both machine-learning applications comes from a Jurassic bimodal carbonate reservoir. Key depositional lithofacies include micritic mudstones, bivalve-coated grain pack/grainstones, Cladocoropsis pack/grainstones, stromatoporoid coralgal pack/grainstones, and skeletal-oolitic grainstones. Key geological features, pore fabric analysis and depositional textures were documented via optical light microscopy. In addition, standard petrophysical measurements (porosity, permeability and grain density) were available from the corresponding core plugs. Although several key uncertainties exist in using machine-learning models for prediction of geological and petrophysical features, these two examples highlight the potential for combining advanced image analysis and data-driven analytical techniques. These approaches can be extended to other multiscale, multidimensional imaging applications in geoscience, such as core and borehole imaging, CT scanning, LIDAR and hyperspectral imaging.@en

The current situation with green gas emission requires the development of low carbon energy solutions. However, a significant part of the modern energy industry still relies on fossil fuels. To combine these two contradictory targets, we investigate a strategy based on a combination of CO2 sequestration with Enhanced Oil Recovery (EOR) in the hydrocarbon reservoirs. In such technology, the development of miscibility is the most attractive strategy from both technological and economic aspects. Modeling of this process involves solving complex nonlinear problem describing compositional flow and transport in highly heterogeneous porous media. An accurate capture of the miscibility development usually requires an extensive number of components to be present in the compositional problem which makes simulation run-time prohibitive for optimization. Here, we apply a multi-scale reconstructing of compositional transport to the optimization of CO2 injection. In this approach, a prolongation operator, based on the parametrization of injection and production tie-lines, is constructed following the fractional flow theory. This operator is tabulated as a function of pressure and pseudocomposition which then is used in the Operator-Based Linearization (OBL) framework for simulation. As a result, a pseudo two-component solution of the multidimensional problem will match the position of trailing and leading shocks of the original problem which helps to accurately predict phase distribution. The reconstructed multicomponent solution can be used then as an effective proxy-model mimicking the behavior of the original multicomponent system. Next, we use this proxy-model in the optimization procedure which helps to improve the performance of the process in several folds. An additional benefit of the proposed methodology is based on the fact that important technological features of CO2 injection process can be captured with lower degrees of freedom which makes the optimization solution more feasible.

@en

As oil production continues worldwide, more oil fields require complex EOR methods to achieve outlined recovery factors. Reservoir engineers are dealing more often with problems involving thermal multiphase multi-component flow models tightly coupled with complex phase behavior. Such modeling implies the solution of governing laws describing mass and energy transfer in the subsurface, which in turn requires the linearization of strongly nonlinear systems of equations. The recently proposed Operator-Based Linearization (OBL) framework suggests an unconventional strategy using the discrete representation of physics. The terms of governing PDEs, discretized in time and space, which depend only on state variables, are approximated by piece-wise multilinear operators. Since the current physical state fully defines operators for a given problem, each operator can be parametrized over the multidimensional space of nonlinear unknowns for a given distribution of supporting points. Onwards, the values of operators, along with their derivatives with respect to nonlinear unknowns, are obtained from the parametrization using multilinear interpolation and are used for Jacobian assembly in the course of a simulation. Previously, the distribution of supporting points was always uniform, requiring higher parametrization resolution to provide accurate and consistent interpolation of an operator around its most nonlinear regions in parameter space. In addition, when the resolution is low, the system can lose hyperbolicity causing convergence issues. In this work, we apply the prior knowledge of underlying physics to distribute the supporting points according to the tie-simplex behavior of the multiphase mixture in parameter space. The approach allows to decrease the parametrization resolution keeping the same accuracy. In addition, the OBL framework is extended to describe multisegment wells working under different controls. We test the accuracy of the developed framework for truly multi-component systems of practical interest.

@en

The nonlinear nature of flow and transport in porous media requires a linearization of the governing numerical-model equations. We propose a new linearization approach and apply it to complex thermal/compositional problems. The key idea of the approach is the transformation of discretized mass- and energy-conservation equations to an operator form with separate space-dependent and statedependent components. The state-dependent operators are parameterized using a uniformly distributed mesh in parameter space. Multilinear interpolation is used during simulation for a continuous reconstruction of state-dependent operators that are used in the assembly of the Jacobian and residual of the nonlinear problem. This approach approximates exact physics of a simulation problem, which is similar to an approximate representation of space and time discretization performed in conventional simulation. Maintaining control of the error in approximate physics, we perform an adaptive parameterization to improve the performance and flexibility of the method. In addition, we extend the method to compositional problems with buoyancy. We demonstrate the robustness and convergence of the approach using problems of practical interest.

@en

Numerical simulation is one of the most important tools required for financial and operational management of geothermal reservoirs. The modern geothermal industry is challenged to run large ensembles of numerical models for uncertainty analysis, causing simulation performance to become a critical issue. Geothermal reservoir modeling requires the solution of governing equations describing the conservation of mass and energy. The robust, accurate and computationally efficient implementation of this solution suggests an implicit time-approximation scheme, which introduces nonlinearity into the system of equations to be solved. The most commonly used approach to solving the system of nonlinear equations is based on Newton's method and involves linearization with respect to nonlinear unknowns. This stage is the most complicated for implementation and usually becomes the source of various errors. A new linearization approach – operator-based linearization – was recently proposed for non-isothermal flow and transport. The governing equations, discretized in space and time, were transformed to the operator form where each term of the equation was specified as the product of two operators. The first operator comprises physical properties of rock and fluids, such as density or viscosity, which depend only on the current state of a grid block, fully defined by the values of nonlinear unknowns. The second operator includes all terms that were not included in the first operators, and depends on both the state and spatial position of a control volume. Next, the first type of operators was parametrized over the physical space of a simulation problem. The representation of highly nonlinear physics was achieved by using multi-linear interpolation, which replaces the continuous representation of parametrized operators. The linearization of the second type of operators was applied in the conventional manner. In this work, we investigated the applicability of this approach to the geothermal processes, specifically for low-enthalpy and high-enthalpy geothermal doublet models with hydrocarbon co-production. The performance and robustness of the new method were tested against the conventional approach on a geothermal reservoir of practical interest. This approach shows significant improvement of geothermal simulation performance, while errors, introduced by coarsening in physics, remain under control. The simplicity of implementation on emerging computational architectures and nonlinearity reduction provide advanced opportunities for uncertainty quantification and risk analysis of geothermal projects.

@en

Reservoir simulation plays an important role in the petroleum industry. Today, there is a specific demand to run ensembles of megaand even giga-cell models. The iterative solution of the large systems of nonlinear governing equations, which describe the multiphase mass transfer in the subsurface, takes the most of the simulation time. The linearization part of the solution process occupies a significant fraction of that time, especially in compositional models. Moreover, the implementation of the linearization step usually embodies the most substantial, complex, and specific part of the computational loop in modern simulators, defining which physical mechanisms and assumptions are employed. This significantly complicates the implementation of simulation codes for heterogeneous computing hardware, which promises significant improvements in simulation time. In this work, we use the recently proposed Operator-Based Linearization (OBL) approach to develop a general purpose reservoir simulation code aiming to substantially decrease the simulation time. OBL offers a simplified linearization method, enhancing the computational performance of simulation and providing an opportunity of a painless porting to heterogeneous computing architectures. To distinct the contribution of both factors, we developed two versions of the compositional simulation prototype code: for traditional CPU and GPU-accelerated hardware architectures. While the former allowed us to speed the linearization stage up by an order of magnitude in comparison with the conventional approach, based on Automatic Differentiation (AD), the latter improved it further by another order of magnitude. The developed prototype realizes the potential of the OBL approach and GPU computing architecture, proving significant improvement in general purpose simulation performance.@en

Abstract The nonlinear nature of flow and transport in porous media requires a linearization of the governing numerical model equations. We propose a new linearization approach and apply it to complex thermal-compositional problems of practical interest. The approach approximates exact physics of a simulation problem, which is similar to an approximate representation of space and time discretization performed in conventional simulation. Maintaining control of the error in approximate physics, we perform an adaptive parametrization to improve the performance of the method and adaptive coarsening to identify the major nonlinearities in the physical description. The global under-relaxation based on inflection points detection in parameterized space of the state variables provides a highly robust strategy to obtain nonlinear solution. This strategy was applied to several challenging compositional problems. In all cases, the proposed strategy improves the robustness and efficiency of the nonlinear solution.@en

General purpose reservoir simulation is based on the solution of governing equations describing mass and energy transfer in the subsurface. The solution process requires the linearization of strongly nonlinear governing equations. Usually, a Newton-based method is used for the linearization. This method demands the assembly of Jacobian and residual for a fully coupled system of equations. Recently, a new linearization approach was proposed and tested for binary systems. The key idea of the Operator Based Linearization (OBL) approach is to transform the discretized mass and energy conservation equations to an operator form which separates space-dependent and state-dependent properties of governing equations. This transformation provides the opportunity to approximate the representation of exact physics (physical properties) of a problem. Specifically, each term of the conservation equations is presented as the product of two different operators. The first operator depends on the current physical state of a system and contains fluid properties, such as density, viscosity, relative permeability, etc. The second operator captures both spatially altered properties, such as permeability, and the rest of state variables, such as pressure in the discrete approximation of gradient. All state-dependent operators are uniformly parametrized within the physical space of the problem (pressure-composition intervals). During simulation process, a multi-linear interpolation is applied to approximate the first type of operators, while the second type of operators is processed based on conventional approach. In this work, we extended the approach to thermal systems with an arbitrary number of components. Besides, we significantly improved the performance of OBL employing adaptive parametrization technique. We tested the approach for truly multi-component thermal systems of practical interest. The computational performance, accuracy, and robustness of a new method were demonstrated against the conventional approach.

@en

Simulation of geothermal processes is based on the solution of strongly nonlinear governing equations describing flow of mass and energy in the subsurface. The solution of this problem requires a linearization of governing equations. Recently, a new approach based on the operator-based multi-linear representation of conservation equations was proposed for an isothermal flow and transport. The discretized version of conservation equations was transformed to the operator form where each term of the equation was represented as a product of two operators: one dependent on physical properties of fluids (e.g. density or viscosity) and another dependent on properties, altered in space (e.g. porosity or permeability). The first type of operators was parametrized in the physical space of the simulation problem at the pre-processing stage, or adaptively during the course of simulation. Next, multi-linear interpolation was applied for physics-based operators when the linearization of the second type of operators was performed in a conventional manner. In this work, the approach has been extended for the thermal processes. The operator-based parametrization was performed using temperature, pressure and composition. The performance and robustness of the method was tested against the conventional approach on low enthalpy models of geothermal reservoirs of practical interest. @en

Non-isothermal multiphase compositional simulation is based on the solution of governing equations describing mass and energy transfer in the subsurface. The solution strategy requires a linearization of strongly nonlinear governing equations describing the process. Usually, a Newton-based method is used for the linearization that demands an assembly of a Jacobian matrix and residuals for a fully coupled system of equations. Recently, a new linearization approach was proposed for compositional problems and tested for simulation of binary compositional and low-enthalpy geothermal flow. The key idea of the approach is the transformation of discretised mass conservation equations to an operator form with separate space-dependent and state-dependent components. This transformation provides an opportunity for an approximate representation of exact physics (physical properties) of the problem. Specifically, each term of conservation equations is represented as a product of two different operators. The first operator depends on a current physical state of a system and contains different properties such as density, viscosity, relative permeability, etc. The second operator captures both spatially altered properties such as permeability and the rest of state variables such as pressure in the discrete approximation of the gradient. At the pre-processing stage, all state-dependent operators are uniformly parametrized within the physical space of the problem (pressure-composition intervals). During the simulation process, a multi-linear interpolation is applied to approximate the first type of operators, while the second type of operators is processed based on the conventional approach. In this work, we have extended this approach to general purpose simulation. We introduced the operator-based parametrization of mass and energy conservations equation based on the pressure, composition, temperature, and porosity. In addition, the approach has been extended and tested on truly multi-component systems of practical interest. The accuracy and robustness of the new method have been tested against the results of simulations based on the conventional approach.@en