Physics-informed Machine Learning for Real-time Unconventional Res- ervoir Management Maruti K. Mudunuru, Daniel O’Malley, Shriram Srinivasan, Jeffrey D. Hyman, Matthew R. Sweeney, Luke Frash, Bill Carey, Michael R. Gross, Nathan J. Welch, Satish Karra, Velimir V. Vesselinov, Qinjun Kang, Hongwu Xu, Rajesh J. Pawar, Tim Carr, Liwei Li, George D. Guthrie, Hari S. Viswanathan Earth and Environmental Sciences Division, Los Alamos National Laboratory, Los Alamos, NM-87544. Department of Geology & Geography, West Virginia University, Morgantown, WV 26506. {maruti,omalled,shrirams,jhyman,msweeney2796,lfrash,bcarey,michael_gross,nwelch,satkarra,vvv,qkang,hxu,ra- jesh,geo,viswana}@lanl.gov, {tim.carr,liwei.li}@mail.wvu.edu Abstract rock formations that have high porosity and permeability We present a physics-informed machine learning (PIML) (Bahadori, 2017). These rock formations are found below an workflow for real-time unconventional reservoir manage- impermeable rock. However, energy extraction from uncon- ment. Reduced-order physics and high-fidelity physics model ventional hydrocarbon resources (Ahmmed and Meehan, simulations, lab-scale and sparse field-scale data, and ma- 2016) involves using advanced drilling and stimulation chine learning (ML) models are developed and combined for real-time forecasting through this PIML workflow. These techniques (e.g., long horizontal laterals and multi-stage hy- forecasts include total cumulative production (e.g., gas, wa- draulic fracturing) to extract crude oil and natural gas that ter), production rate, stage-specific production, and spatial are trapped in the pores of relatively impermeable sediment- evolution of quantities of interest (e.g., residual gas, reservoir ary rocks (e.g., shale, tight sandstones). pressure, temperature, stress fields). The proposed PIML Typically, unconventional reservoirs have porosity in the workflow consists of three key ingredients: (1) site behavior libraries based on fast and accurate physics, (2) ML-based range of 0.04-0.08 and matrix permeability on the order of inverse models to refine key site parameters, and (3) a fast nanodarcies (10-16-10-20 m2) (Rezaee, 2015; Belyadi et al., forward model that combines physical models and ML to 2019). Instead of the porous flow that dominates conven- forecast production and reservoir conditions. First, synthetic tional reservoirs, fracture flow dominates unconventional production data from multi-fidelity physics models are inte- reservoirs, with natural fractures dissecting the matrix and grated to develop the site behavior library. Second, ML-based inverse models are developed to infer site conditions and en- intersecting with the hydraulic fractures. As result, energy able the forecasting of production behavior. Our preliminary extraction is more difficult than conventional reservoirs. results show that the ML-models developed based on PIML Model-based optimization of unconventional reservoirs is workflow have good quantitative predictions (>90% based on also challenging because due to the long horizontal laterals R2-score). In terms of computational cost, the proposed ML- there is insufficient site data to inform high-fidelity physics models are π’ͺπ’ͺ(104 ) to π’ͺπ’ͺ(107 ) times faster than running a high-fidelity physics model simulation for evaluating the models (Mohaghegn, 2017; Belyadi et al., 2019). Despite quantities of interest (e.g., gas production). This low compu- these challenges and due to the abundance of unconven- tational cost makes the proposed ML-models attractive for tional resources, with reserves projected to last for many real-time history matching and forecasting at shale-gas sites decades, energy extraction from these resources have gained (e.g., MSEEL – Marcellus Shale Energy and Environmental prominence in recent years (Briefing, 2013 and Weijermars, Laboratory) as they are significantly faster yet provide accu- rate predictions. 2014). Current extraction efficiency from unconventional reservoir is very low (~5-10%) for tight oil and ~20% for 1. Introduction shale gas (Sandrea, 2007; Muggeridge et al., 2014) com- Energy extraction from conventional resources involves pared to conventional reservoirs (~20-40%) (Zitha et al., producing crude oil, natural gas, and its condensates from 2008). This is because the impact of resource development Copyright Β© 2020, for this paper by its authors. Use permitted under Crea- tive Commons License Attribution 4.0 International (CCBY 4.0). processes (e.g., slow drawdown or fast drawdown) and un- training regimes, and the exploration of novel production derlying physics that determine the energy extraction from strategies fundamentally requires extrapolation (where ML the impervious rocks are poorly understood (Rezaee, 2015; struggles) as opposed to interpolation (where ML excels). Belyadi et. Al., 2019). State-of-the-art workflows for uncon- The physics-based workflows adopted for modeling con- ventional reservoir management are data-driven, which per- ventional reservoirs use extensively available site-character- form poorly beyond their training regimes. Hence, innova- ization data (which is acquired over months or years) for tive extraction strategies (e.g., pressure-drawdown manage- history-matching. Even though large unconventional reser- ment) coupled with advanced workflows (e.g., physics-in- voir data (e.g, fiber-optics) are sampled at sparse locations, formed machine learning) are needed to improve the hydro- simply combining all the data together would not improve carbon recovery efficiency (Seales et al., 2017; Lougheed et the accuracy of real-time forecasting. This is because reser- al., 2017; Mirani et al., 2018). In this paper, we present a voir conditions change considerably from one basin to an- physics-informed machine learning (PIML) workflow other basin. Therefore, physical constraints need to be incor- (Fig.1) to address unconventional production for real-time porated in workflows. Existing workflows employ high-fi- reservoir management. One of the goals of this PIML work- delity physics models to perform simulations, which are ex- flow (Fig.2) is to develop fast and accurate ML-models pensive to run. For example, it takes several days to months grounded in physics for real-time history matching and pro- to run reservoir-scale model simulations (Rezaee, 2015; duction forecasting in a fracture shale gas reservoir. Rajput and Thakur, 2016; Belyadi et al., 2019) with degrees- of-freedom in the π’ͺπ’ͺ(108 ) on state-of-the-art HPC machines. 1.1 State-of-the-art Workflows and Key Gaps As a result, these conventional reservoir workflows are not Current workflows for unconventional reservoirs are pre- ideal for usage in comprehensive uncertainty quantification dominantly based on production decline curve analysis and studies, which require 1000s of forward model runs. To its extensions (Wu et al., 2013; Sun 2015), data-driven ma- overcome the key gaps associated with existing workflows, chine learning (ML) approaches (Holdaway 2014; Chaki we propose a PIML workflow (Fig.1 and Fig.2) to accelerate 2015; Puyang et al., 2015; Carvajal et al. 2017; Mohaghegn, the development of ML-models while constraining it with 2017), and/or extension of physics-based conventional res- physics. The aim is (1) to develop a library from a combina- ervoir workflows (Rezaee, 2015; Rajput and Thakur, 2016; tion of site observations and physics-based models that is Belyadi et al., 2019). Decline curve analysis provides em- representative of unconventional reservoir site behavior, pirical models to forecast production data based on the past and (2) to develop fast and accurate ML-models for real- production history. This type of approach lacks physics and time history matching to refine key site parameters and fore- in-depth knowledge of the site behavior is not included in cast production quantities of interest (QoIs) with uncertainty the forecasting models. The data-driven ML approaches per- estimates. Key production QoIs include, total cumulative form poorly when faced with uncertain, missing, and sparse production of hydrocarbons, gas and water production rates, data – common problems with existing datasets related to stage-specific production, and spatial evolution of residual unconventional reservoirs. Moreover, the data-driven ML gas, reservoir pressure, temperature, and stress fields based analyses perform poorly in making forecasts outside of their on a user-defined pressure-drawdown strategy. Figure-1: Physics-informed machine learning (PIML) workflow for reservoir management. inverse model, then the smaller number of runs from the high-fidelity model are used to fine-tune the inverse model. This effectively represents a multi-fidelity ap- 2. Physics-informed Machine Learning proach to training the ML inverse model. This ML-in- 2.1 Innovation and proposed approach verse model provides capabilities for real-time history matching to update the key parameters (e.g., rock perme- The proposed approach to develop the PIML workflow ability, rock porosity, gas transport properties). Sensitiv- consists of three key steps: (Step-1) Development of site ity analysis (e.g., Sobel indices, Random Forests) is per- behavior libraries based on fast and accurate physics, formed to provide quantitative information on the key (Step-2) Development of ML-based inverse models to in- sensitive parameters (e.g, matrix permeability and poros- fer key site parameters, and (Step-3) Development of fast ity, fracture network parameters) that influence shale gas forward models that combine physical models with ML production rates. to forecast production and reservoir conditions. PIML Workflow Step-3: The physics-based reduced PIML Workflow Step-1: Development of a site be- order forward model uses these calibrated parameters for havior library involves generating synthetic data for a real-time forecasting. This model is trained on the site li- range of possible site characteristics. This includes a large braries along with evolving production data. Moreover, it number of runs from a fast physics-based reduced-order allows for various operational decisions (e.g, slow draw- model and a smaller number of runs from a high-fidelity, down vs. fast drawdown) to be evaluated relative to future full physics model. The fast physics models (e.g., Patzek outcomes. The ML-inverse and physics-based reduced models) allow us to quickly model and build a site library order forward model can be combined to provide uncer- on the evolving reservoir data. These models allow us to tainty estimates on the production quantities of interest efficiently explore the parameter space. Moreover, com- (e.g., remaining hydrocarbon-in-place, spatial evolution bining fast physics models with ML allows us to identify of pressure and temperature, shale gas and water produc- the important physical processes or dominant mecha- tion as a function of time). nisms that must be represented during full physics simu- lations. The dominant mechanisms at different stages of 2.2 Key Challenges production include The key challenges to accelerate the proposed PIML work- β€’ Early stages of production (<1 year): Primary frac- flow include: ture creation, geometry, and network connectivity, 1Q. How do we enhance the information content and fill the primary fracture behavior, propped fracture behavior, gaps in the limited unconventional site-data for PIML anisotropic permeability of fractures. analyses? Specifically, how to use the short-time gas pro- β€’ Middle stages of production (~1-5 years): Second- duction data (e.g., 30-120 days) to forecast the long-term ary fractures and their permeabilities, shear fracture performance (e.g., 1-5 years) of an unconventional shale geometries, and geochemical impacts of hydraulic flu- gas well? ids and formation water. 2Q. What is the minimum number of high-fidelity physics β€’ Late stages of production (~5-10 years): Matrix po- model simulations (e.g., PFLOTRAN, FEHM, dfnWorks) rosity and transport properties, water imbibition im- needed to develop a gas production library that is repre- pacts, adsorbed gas properties, pore structure distribu- sentative of shale gas sites behavior? How do we improve tion. the full physics models to accurately represent the site be- These mechanisms are not accounted for in the fast phys- havior? ics models but are needed to describe the reservoir behav- 3Q. How can we use the information learned from reduced- ior at different stages of gas and water production. How- order physics models (e.g., Patzek model, graph-based ever, simulating these mechanisms using full physics models) to inform high-fidelity physics model simula- models for entire parameter space is computationally in- tions? tractable. As a result, a fast physics model library along 4Q. What are the key sensitive parameters (e.g., matrix and with ML are used to guide and improve the development fracture properties) in high-fidelity physics models that of the full physics model library (which consist of complex influence short-term and long-term gas production rates? 3D simulations of matrix-fracture interactions) for char- 5Q. What are the key operational parameters (e.g., pressure acterizing site behavior. drawdown strategies) that can be used to inform decisions PIML Workflow Step-2: The ML-inverse model is in real-time, leading to optimized production? developed on this site behavior library using a transfer learning approach (Pan et al., 2009; Goodfellow et al., 2.3 Our Hypothesis 2016; Yamada et al., 2018) where the numerous runs Our hypotheses/approaches to address the key challenges from the reduced order model are used to train an initial include: 1A. Augment limited unconventional reservoir data (e.g., PetroMehras) can be used to generate synthetic data for the short-term production data, well logs) with lab-scale ex- site behavior library. perimental data, reduced-order physics simulation data, The reduced-order physics models (Patzek et al., 2013) to and high-fidelity physics simulation data. Short-term pro- generate synthetic data are given by duction data (e.g., 30-120 days) contains statistical infor- πœ•πœ•π‘šπ‘š οΏ½ οΏ½ ) πœ•πœ•π‘šπ‘š πœ•πœ• 𝛼𝛼(π‘šπ‘š οΏ½ mation (e.g., matrix and fracture network properties) that βˆ’ οΏ½ οΏ½=0 (1) Μƒ πœ•πœ•π‘‘π‘‘ πœ•πœ•π‘₯π‘₯οΏ½ 𝛼𝛼𝑖𝑖 πœ•πœ•π‘₯π‘₯οΏ½ can help us predict the long-term production data (e.g., 1- 5 years). where 𝑑𝑑̃ is the dimensionless time, π‘₯π‘₯οΏ½ is the dimensionless 2A. Reduced-order physics models provide information on distance, and π‘šπ‘š οΏ½ is the real gas pseudopressure. Eq.(1) cor- the key sensitive parameters needed to inform high-fidel- responds to reduced-order physics of gas flow in fractured ity physics model simulations. As we obtain more produc- porous rock. These are given as follows tion data and/or site data, a ML-based active-feedback 𝑑𝑑 𝑑𝑑 2 π‘₯π‘₯ π‘˜π‘˜ loop is used to improve full physics models. In this active 𝑑𝑑̃ = 𝜏𝜏 = π‘₯π‘₯οΏ½ = 𝛼𝛼𝑖𝑖 = οΏ½ 𝜏𝜏 𝛼𝛼𝑖𝑖 𝑑𝑑 πœ™πœ™π‘ π‘ π‘”π‘” πœ‡πœ‡π‘”π‘” 𝑐𝑐𝑔𝑔 𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 𝑝𝑝,𝑇𝑇 feedback loop, ML-inverse model is used to update and constrain the key parameters based on newly available 1 �𝑐𝑐𝑔𝑔 π‘π‘οΏ½πœ‡πœ‡π‘”π‘” 𝑍𝑍𝑔𝑔 reservoir data. Based on these updated key parameters, π‘šπ‘š οΏ½= οΏ½ οΏ½ π‘šπ‘š(π‘₯π‘₯, 𝑑𝑑) 2 𝑝𝑝2 site-behavior libraries are also updated to account for new 𝑖𝑖 production data. 𝑝𝑝 3A. Transfer learning can be used to provide link and trans- 𝑝𝑝 𝑑𝑑𝑑𝑑 π‘šπ‘š(𝑝𝑝, π‘₯π‘₯, 𝑑𝑑) = 2 οΏ½ fer information from reduced-order physics to high-fidel- πœ‡πœ‡π‘”π‘” 𝑍𝑍𝑔𝑔 (𝑝𝑝, 𝑇𝑇) ity physics. π‘π‘π‘π‘β„Žπ‘π‘ 4A. Sensitivity analysis (e.g., Sobel indices, Random For- ests) can provide quantitative information on the key sen- π‘˜π‘˜(𝑝𝑝) (2) π›Όπ›ΌοΏ½π‘šπ‘š(𝑝𝑝)οΏ½ = οΏ½ sitive parameters (e.g, matrix permeability, matrix lithol- οΏ½πœ™πœ™π‘ π‘ π‘”π‘” + (1 βˆ’ πœ™πœ™)π’¦π’¦π‘Žπ‘Ž οΏ½πœ‡πœ‡π‘”π‘” 𝑐𝑐𝑔𝑔 𝐸𝐸𝐸𝐸𝐸𝐸𝐸𝐸𝐸𝐸𝐸𝐸𝐸𝐸𝐸𝐸 𝑝𝑝,𝑇𝑇 ogy, fracture length and orientation, stage spacing, hydro- carbon-in-place) that influence short-term and long-term where 𝜏𝜏 is the characteristic inference time, 𝑑𝑑 is the half- gas production rates through feature importance. distance between the hydrofractures, 𝛼𝛼𝑖𝑖 is the initial hydrau- 5A. Maximizing early production may not maximize total lic diffusivity, π‘˜π‘˜ is the fractured rock effective permeability recovery efficiency. Optimal pressure management strat- dependent on reservoir pressure, πœ™πœ™ is the fractured rock ef- egies are needed to enhance recovery efficiency. One of fective porosity, 𝑠𝑠𝑔𝑔 is the fraction of pore space occupied by our hypotheses is that slower drawdown rates can lead to the gas, πœ‡πœ‡π‘”π‘” is the gas viscosity, 𝑐𝑐𝑔𝑔 is the isothermal com- improved recovery efficiency in long-term shale gas pro- pressibility of gas, π’¦π’¦π‘Žπ‘Ž is the differential equilibrium parti- duction (e.g., 5-10 years). tioning coefficient of gas, and 𝑍𝑍𝑔𝑔 is the compressibility fac- tor of gas. These gas properties are dependent on evolving 2.4 Details of the Proposed Approach reservoir pressure and temperature. Eq. (1) is a nonlinear Fig.1 shows the overall PIML workflow for reservoir pressure diffusion equation, which is solved numerically. management. Short-time production data is fed to ML-in- The cumulative production of gas mass is given as follows verse model to perform history-matching and infer key site 𝑑𝑑̃ parameters. These key parameters are then fed to ML-for- πœ•πœ•π‘šπ‘š οΏ½ 𝔐𝔐(𝑑𝑑̃) = β„³ οΏ½ οΏ½ 𝑑𝑑𝑑𝑑 β€² (3) ward model to forecast long-term production QoIs. Fig.2 πœ•πœ•π‘₯π‘₯οΏ½ π‘₯π‘₯οΏ½=0 0 and Fig.3 provide more details on our PIML workflow. where β„³ is the initial hydrocarbon-in-place. These figures show how the site behavior library is con- The expensive full physics models to simulate gas flow structed and used to develop ML-models for history match- and transport in fractured porous media (Rezaee, 2015; ing and real-time forecasting. Physically-realistic synthetic Salama et al., 2017; Belyadi et al., 2019) are given by data is generated for a range of possible site characteristics πœ•πœ•πœ•πœ•π‘ π‘ π‘”π‘” πœŒπœŒπ‘”π‘” using multi-fidelity physics models. This synthetic data pro- + 𝛻𝛻 βˆ™ πœŒπœŒπ‘”π‘” 𝒒𝒒 = 0 (4) vides insights on relevant features on poorly constrained site πœ•πœ•πœ•πœ• parameters and production scenarios that are not-yet ob- π‘˜π‘˜(𝑝𝑝) 𝒒𝒒 = βˆ’ 𝛻𝛻�𝑝𝑝 βˆ’ πœŒπœŒπ‘”π‘” 𝑔𝑔𝑔𝑔� (5) served. The site behavior library can be updated actively as πœ‡πœ‡π‘”π‘” the field-scale data becomes available overtime. This allows πœ•πœ•πœ•πœ•πœ•πœ• (6) + 𝛻𝛻 βˆ™ (𝑐𝑐𝒒𝒒 βˆ’ πœ™πœ™πœ™πœ™πœ™πœ™πœ™πœ™πœ™πœ™πœ™πœ™) = 0 us to update or prune parameter combinations that are in- πœ•πœ•πœ•πœ• consistent with site characteristics. Note that any simulation where πœŒπœŒπ‘”π‘” is the density of the gas which is dependent on the platform (e.g., ECLIPSE, INTERSECT, tNavigator, CMG reservoir pressure, 𝒒𝒒 is the Darcy flux, 𝑐𝑐 is the gas concen- suite, Landmark Nexus, MRST, BOAST, OPM) (see ref. tration, πœ†πœ† is the tortuosity, and 𝐷𝐷 is the fracture rock effec- model. The loss function for this ML-inverse model is de- tive diffusivity. Eq.(4) and Eq.(5) model the gas flow under fined in terms of how well it works in combination with the varying reservoir and bottom hole pressures. The underlying reduced-order physics model at predicting future produc- assumptions include Darcy’s flow and Fick’s law. Adsorp- tion. Note that the second library can also be augmented tion and non-pore refinement effects on phase behavior are with field production data to improve the realism of the ML- ignored. Eq.(6) models the gas transport from fractured inverse model. This fine-tuning represents a multi-fidelity stage to horizontal well based on the initial hydrocarbon-in- approach to machine learning where a large dataset is gen- place. The amount of gas extracted from the pair of hydro- erated with a reduced-order physics model and a smaller da- fractures at a given section of the well is equal to 𝑐𝑐𝒒𝒒. Differ- taset is generated with a high-fidelity, expensive full physics ent types of equation of state (EOS) models are used to eval- model. This multi-fidelity approach allows us to perform uate the gas density. These include ideal gas law, exponen- real-time history matching on new production data to deter- tial pressure-dependent model, and Redlich-Kwong-Soave mine critical site parameters that can be used to accurately model. The corresponding EOS model expressions are given predict future production. by 𝑝𝑝𝑀𝑀𝑔𝑔 πœŒπœŒπ‘”π‘” = 3. Results 𝑅𝑅𝑅𝑅 Fig.4 shows a DFN model of a single stage based on field πœŒπœŒπ‘”π‘” = 𝜌𝜌0 𝑒𝑒 𝛽𝛽×(π‘π‘βˆ’π‘π‘0) data from the Marcellus Shale Energy and Environment La- 𝑝𝑝𝑀𝑀𝑔𝑔 boratory (MSEEL) shale gas site. This model was built us- πœŒπœŒπ‘”π‘” = (7) ing dfnWorks, which is a computationally expensive full 𝑍𝑍𝑔𝑔 (𝑝𝑝, 𝑇𝑇)𝑅𝑅𝑅𝑅 physics software suite used to generate high-resolution rep- Note that it is possible to incorporate nanopore confinement resentations of DFNs (Hyman et al., 2015). This high-fidel- effects (e.g,shifts in bubble or dew points) into EOS models ity meshing of the fracture network is critical to accurately to account for density changes. For example, see ref. Islam capture the physical processes in a fractured shale gas reser- et al., 2015, Tan and Piri, 2015, and Liu and Zhang, 2019. voir. To capture the matrix effects, we generate an octree- Through these multi-fidelity physics models, the reser- refined continuum mesh (grey color in Fig.4) based on the voir physical behavior is captured accurately. Gas flow and DFN. The original DFN model in Fig.4 consists of three hy- transport mechanisms are accounted through conservation draulic fractures and a swarm of natural fractures that are of mass and equation of state for real gases. State-of-the-art connected to hydraulic fractures. While generating the con- simulators (e.g., PFLOTRAN, dfnWorks) are used to de- tinuum mesh, the DFN model is simultaneously upscaled to velop high-fidelity simulation data. These simulators use fi- account for matrix-fracture interactions, which results in ac- nite volume methods and Newton-Raphson method to solve curate permeability and porosity values that are needed to the discretized system of nonlinear equations given by simulate gas flow and transport in fractured shale. The final Eq.(4)-(7). Moreover, these simulators account for accurate mesh contains approximately 500,000 mesh cells. meshing of fractures, matrix, and upscaling of fracture net- Fig.5 shows the flow and transport simulation using work properties for reservoir-scale high-fidelity physics PFLOTRAN with a Barton-Bandis stress relationship (Bar- simulations. ton and Bandis, 1990). This figure shows the drainage over Fig.3 shows the PIML workflow to create efficient in- a period of 10 years. For gas flow simulations (left), the well verse models to infer key site parameters from production is maintained at 12MPa and the reservoir initial pressure is information. First, we sample the relevant regions in param- at 21MPa. For gas transport simulations (right), the initial eter space. Two site behavior libraries are developed based hydrocarbon-in-place is assumed to spread in the entire frac- on this sampling. One comes from running the reduced-or- ture stage and the figure shows the transport of hydrocarbon der physics model with a large set of samples and the other to the well over time at two different vertical heights. The comes from running the full physics model with a smaller main inference from these simulations is that the character- set of samples. The first library (based on the reduced-order istic behavior of drainage is tied to hydraulic fractures. Our physics model) is used to train an initial ML-inverse model. future work involves accounting for heterogeneity of hydro- The ML-inverse model is then fine-tuned with the library carbon distribution in the matrix for production forecasts. from the high-fidelity physics model using transfer learning, Fig.6 shows encouraging results of long-term production producing a final ML-inverse model. During this fine-tuning forecasts using ML-forward models. The red color repre- process the weights of the neural network are fine-tuned. sents the short-term production, which is used along with The ML-inverse model takes past production as input and the site behavior library built on reduced-order physics mod- produces physical parameters as output. These physical pa- els to infer key site parameters (e.g., initial hydrocarbon-in- rameters can then be fed into the reduced-order physics place, hydraulic diffusivity). History-matching is performed on the short-term production data. That is, the ML-inverse full physics models. The initial ML-inverse model trained model is trained on 90 days of production data (red color). on the reduced physics site library and short-term produc- Then, the resulting ML-model with site behavior library tion data provided us key site parameters, which are hydrau- based on reduced-order physics is used to predict future pro- lic diffusivity and initial hydrocarbon-in-place. This initial duction (green color). The markers represent the long-term ML-inverse model is then fine-tuned with the library from production data and the solid green color line represents the the expensive full physics models using transfer learning. ML-model forecasts. Quantitatively, the prediction accu- The expensive full physics model simulations were devel- racy based on R2-score is > 90%. From this figure, it is clear oped using dfnWorks and PFLOTRAN simulators. These that ML-model predictions, which are combined with phys- high-fidelity simulations account for matrix-fracture inter- ics are able to accurately represent the long-term production action, which is needed to accurately simulate gas flow and data. transport in fractured shale reservoirs. Moreover, from these simulations we inferred that the characteristic behavior of 3.1 Discussion drainage is tied to hydraulic fractures. This high-fidelity We note that if another site/formation (e.g., Woodford, simulation library was used to fine-tune the ML-inverse Haynesville, Fayetteville, Barnett, Utica, EagleFord) has a model. Our ongoing work seeks to advance and test the similar parameter range, we can use a set of ML techniques PIML workflow, site behavior libraries, and ML-models for derived from transfer learning to model this site/formation. pressure management (e.g., slow drawdown vs. fast draw- Transfer learning helps us to transfer knowledge gained down) to optimize recovery at MSEEL. Our preliminary re- from one site (e.g., Marcellus) to another site (e.g., Wood- sults shown in this paper are encouraging for use of site be- ford, Barnett). But the developed ML-models for one site havior libraries with ML-inverse model to address this prob- (e.g., MSEEL) need fine-tuning (or minimal retraining the lem. neural networks) to transfer knowledge across shale sites/formations. This is not burdensome when compared to developing a new site behavior library and ML-model for a Acknowledgements different site altogether. The transfer learning approach is This work was funded by the U.S. Department of Energy attractive for tasks where reusability of ML-models for sim- (DOE) – Office of Fossil Energy (FE) through its Oil and ilar types is of great importance. Natural Gas Research Storage Program as implemented by If the production curve contains high frequency content the National Energy Technology Laboratory. We would or oscillations, which might need addressing, there are vari- specifically like to thank Jared Ciferno for his support and ous ways to incorporate this in the ML analyses and physics- guidance. LANL is operated by Triad National Security, based models. For example, from the production curve and LLC, for the National Nuclear Security Administration of bottom hole pressure, we can extract dominant frequencies U.S. Department of Energy (Contract No. or a range of frequencies we are interested in modeling 89233218CNA000001). The opinions expressed in this pa- through Fast Fourier Transformation (FFT). This FFT trans- per are those of the authors and do not necessarily reflect formation of bottom hole pressure and production curve vs. that of the sponsors. time provides in-depth information on quantitative aspects of high-frequency oscillations to be incorporated in physics models (e.g., π‘π‘π‘π‘β„Žπ‘π‘ = π‘Žπ‘Ž1 sin πœ”πœ”1 𝑑𝑑 + π‘Žπ‘Ž2 sin πœ”πœ”2 𝑑𝑑 + References … + π‘Žπ‘Žπ‘›π‘›βˆ’1 sin πœ”πœ”π‘›π‘›βˆ’1 𝑑𝑑 + π‘Žπ‘Žπ‘›π‘› sin πœ”πœ”π‘›π‘› 𝑑𝑑 ) when developing the Ahmmed, U. and Meehan, D. N., (2016), Unconventional Oil and site behavior library and pressure management strategies Gas Resources: Exploitation and Development, CRC Press, Boca (e.g., drawdown frequencies). Raton, FL, USA. Bahadori, A. (2017), Fluid Phase Behavior for Conventional and Unconventional Oil and Gas Reservoirs, Elsevier, Oxford, UK. 4. Conclusions Barton, N., and Bandis, S. (1990), Review of predictive capabili- In this paper, we have presented a PIML workflow for real- ties of JRC-JCS model in engineering practice. In Rock Joints, time history matching and forecasting of gas production Proc. Int. Symp on Rock Joints, Loen, Norway (eds N. Barton and QoIs. The workflow coupled the strengths of machine learn- O. Stephenson) (pp. 603-610). ing with the predictability of physics-based models for real- Belyadi, H., Fathi, E., and Belyadi, F., (2019), Hydraulic Fractur- time history-matching and forecasting. The PIML workflow ing in Unconventional Reservoirs: Theories, Operations, and Eco- nomic Analysis, Gulf Professional Publishing, Elsevier, Cam- used short-term production data and site behavior libraries bridge, Massachusetts, USA. to perform real-time history matching. Site behavior librar- Briefing, U. S. S. (2013), International energy outlook 2013, US ies are developed based on many runs of a reduced-order Energy Information Administration. physics models and a smaller number of runs of expensive Carvajal, G., Maucec, M., and Cullick, S. (2017), Intelligent Digi- Patzek, T. W., Male, F., and Marder, M. (2013), Gas production in tal Oil and Gas Fields: Concepts, Collaboration, and Right-time the Barnett Shale obeys a simple scaling theory. Proceedings of the Decisions. Gulf Professional Publishing, Elsevier, Cambridge, National Academy of Sciences, 110(49), 19731-19736. Massachusetts, USA. PetroMehras -- List of Petroleum Softwares and Simulators, Chaki, S., (2015), Reservoir Characterization: A Machine Learn- https://www.petromehras.com/petroleum-software-directory/res- ing Approach. arXiv preprint arXiv:1506.05070. ervoir-simulation-software/list-of-petroleum-software-application Goodfellow, I., Bengio, Y., and Courville, A. (2016), Deep learn- PFLOTRAN: https://www.pflotran.org/ ing. MIT press. Puyang, P., Taleghani, A. D., and Sarker, B. R., (2015), Multi-dis- Holdaway, K. R. (2014), Harness Oil and Gas Big Data with An- ciplinary data integration for inverse hydraulic fracturing analysis: alytics: Optimize Exploration and Production with Data-driven A case study. SPE Journal, 2015. Models. John Wiley & Sons, New Jeresey, USA. Rajput, S., and Thakur, N. K., (2016), Geological Controls for Gas Hyman, J. D., Karra, S., Makedonska, N., Gable, C. W., Painter, S. Hydrates and Unconventionals, Elsevier, Cambridge, Massachu- L., and Viswanathan, H. S. (2015), dfnWorks: A discrete fracture setts, USA. network framework for modeling subsurface flow and transport, Rezaee, R. (2015), Fundamentals of Shale Gas Reservoirs, John Computers & Geosciences, 84, 10–19. Wiley & Sons, New jersey, 2015. Islam, A. W., Patzek, T. W., and Sun, A. Y. (2015), Thermody- Salama, A., Amin, M. F. E., Kumar, K., and Sun, S. (2017), Flow namics phase changes of nanopore fluids, Journal of Natural Gas and transport in tight and shale formations: A review. Geofluids. Science and Engineering, 25, 134-139. Sandrea, I. & Sandrea, R. (2007), Global oil reserves - 1: Recovery Liu, X., and Zhang, D. (2019), A review of phase behavior simu- factors leave vast target for EOR technologies. Oil Gas Journal. lation of hydrocarbons in confined space: Implications for shale oil and shale gas, Journal of Natural Gas Science and Engineering, Seales, M. B., Ertekin, T., and Wang, J. Y. (2017), Recovery effi- 68, 102901. ciency in hydraulically fractured shale gas reservoirs. Journal of Energy Resources Technology, 139(4), 042901. Lougheed, D., and Anderson, D. (2017), Does Pressure Matter? A Statistical Study. In SPE’s Unconventional Resources Technology Sun, H. (2015), Advanced Production Decline Analysis and Appli- Conference (UreTEC), Austin, Texas (pp. 1919-1933). cation. Gulf Professional Publishing, Elsevier, Cambridge, Massa- chusetts, USA. Mirani, A., Marongiu-Porcu, M., Wang, H., and Enkababian, P. (2018), Production-pressure-drawdown management for fractured Tan, S. P., and Piri, M. (2015), Equation-of-state modeling of con- horizontal wells in shale-gas formations. SPE Reservoir Evalua- fined-fluid phase equilibria in nanopores, Fluid Phase Equilibria, tion & Engineering, 21(03), 550-565. 393, 48-63. Mohaghegh, S. D., (2017), Shale Analytics: Data-driven Analytics Weijermars, R. (2014), U. S. shale gas production outlook based in Unconventional Reservoirs, Springer, Switzerland. on well roll-out rate scenarios, Applied Energy, 124, 283-297. MSEEL: http://mseel.org/ Wu, X., He, J., Zhang, H., & Ling, K. (2013), Tactics and Pitfalls in Production Decline Curve Analysis. In SPE Production and Op- Muggeridge, A., Cockin, A., Webb, K., Frampton, H., Collins, I., erations Symposium. Moulds, T., and Salino, P. (2014), Recovery rates, enhanced oil recovery and technological limits. Philosophical Transactions of Yamada, M., Chen, J., and Chang, Y. (2018), Transfer Learning: the Royal Society A: Mathematical, Physical and Engineering Sci- Algorithms and Applications. Morgan Kaufmann. ences, 372(2006), 20120320. Zitha, P., Felder, R., Zornes, D., Brown, K. & Mohanty, K. (2008), Pan, S. J., and Yang, Q. (2009), A survey on transfer learning. Increasing Hydrocarbon Recovery Factors. SPE Journal. IEEE Transactions on knowledge and data engineering, 22(10), 1345-1359. Figure-2: Details of PIML workflow. The workflow utilizes a library of data on site behavior to inform forecasting models. The scale of physical parameters used in the development of site behavior libraries are stage spacing (~100-200m), hydraulic fracture length (~100-150m), a stage may contain 3-4 hydraulic fractures and a swarm of natural fractures that are connected to hydraulic fractures, fracture rock reservoir permeability (~0.1-0.9mD), fractured rock porosity (~0.04-0.08), reservoir pressures (~20- 30MPa), well flowing pressures (~5-15MPa), and gas properties (e.g., density, saturation, viscosity, compressibility). Figure-3: PIML workflow to create an efficient set of ML-inverse models for real-time history-matching. Figure-4: Discrete fracture network (DFN) and upscaled DFN models to simulate the physical behavior of the fractured reservoir. Figure-5: Gas flow and transport in upscaled DFN model using PFLOTRAN simulator with Barton-Bandis stress relationship. Figure-6: Long-term production forecasts using ML-forward model trained on reduced-order physics models.