Fracture Sustainability in Enhanced Geothermal Systems: Experimental and Modeling Constraints

Several studies [1–3] have highlighted the potential of enhanced geothermal systems (EGS) as a source of energy for generating electricity in the United States from 100 GWe up to 5157 GWe. There are a number of critical technical challenges that need to be addressed to make EGS a technically and economically viable option: these include reservoir access (through improved drilling technology), reservoir creation (through improved fracture stimulation methods), and reservoir sustainability. To address these key research topics, the US Department of Energy’s Geothermal Technologies Office has launched the Frontier Observatory for Research in Geothermal Energy (FORGE) initiative to stimulate R&D advances needed to make EGS viable through the use of a dedicated field site. The FORGE Roadmap [4] identifies three critical research areas for EGS: stimulation planning and design, fracture control, and reservoir management. These topics involve the stimulation, control, and sustainability of fractures needed for circulation of a working fluid to extract heat from an EGS reservoir [5].

There are two main approaches to creating fracture permeability for an EGS reservoir: shear stimulation of critically stressed fractures, and hydraulic fracturing to create new fractures [6]. In order to keep these fractures open after stimulation, they either need to be injected with proppant (as is commonly done with hydraulic fracturing in unconventional oil and gas reservoirs), or propped open by asperities on the fracture surface caused by shear offsets along a rough natural fracture surface. Mixed-mode stimulation, where both shear and tensile failure of fractures occur, is another possible option for increasing fracture permeability in an EGS reservoir [7].

Fractures in natural hydrothermal systems often undergo cycles of opening and closing related to stress-induced dilation, mineral precipitation, brecciation, and local dissolution and precipitation [15,16]. Davatzes and Hickman [15] observe that fractures containing ductile sheet silicate minerals, such as smectite and chlorite, are much less likely to remain open compared with those containing more brittle mineral phases such as quartz and calcite. Thus, the nature of both the host rock mineralogy as well as minerals precipitated along the fracture surfaces can play an important role in the relative reactivity of the fracture surface as well as its geomechanical stability—this is especially important for those phases forming the fracture asperities. These fracture surface features play an important role in preferential dissolution, self-propping, and the creation of new asperities [17].

The objective of this study is to evaluate the longevity of fracture apertures (and thus the sustainability of fluid circulation within an EGS reservoir) through a series of hydrothermal experiments using fractured rock samples of varying compositions at geothermally relevant pressure and temperature conditions. These experiments were performed using a purpose-built apparatus that allows the application of a normal stress to fractured core. Core samples were selected from two geothermal fields where infield EGS field tests have been conducted (Desert Peak [18–21] and Bradys [21]) as well as a granite, a rock type relevant to the Utah FORGE EGS field site [22]—these samples provide a range of lithologic types for potential EGS reservoirs. These well-characterized tests were used to constrain THMC numerical models used to simulate one of the experiments. The results of the laboratory experiments and numerical simulations are presented in Secs. 2 and 3.

The hydrothermal-mechanical experiments were conducted in LBNL’s Geothermal Laboratory using a specially designed apparatus to evaluate how normal stress, fluid pressure, and water–rock interaction at geothermally relevant temperatures might change the flow behavior of an asperity-propped fracture. These tests were run using four distinct rock types. Descriptions of the test apparatus, the rock samples selected for these experiments, and the experimental procedures are given below.

A custom-built uniaxial stress frame with a heater and flow system was designed and created for these experiments (Fig. 1). The purpose of this system is to flow water through a fractured rock core at temperatures and normal stresses relevant to EGS reservoirs over time and detect any changes in fluid flow pressure that would reflect a change in fracture permeability. The system allows for uniaxial stress to be applied to a rock core with a fracture that has been created oriented perpendicular to the core axis. To allow water to flow into the fracture, a small diameter hole was drilled parallel to the core axis through where the fracture was created to allow fluid to be introduced into the center of the fracture—the resulting flow would be nominally radial along the fracture face to the edge of the core.

The pressure vessel was constructed of Grade 5 titanium (containing ∼6% Al) that was cleaned and passivated by soaking in concentrated nitric acid and heated in an oven at 400 °C for 12 h to create an oxidized layer to reduce potential corrosion and contamination of the injected fluid. The system was sealed at the top and bottom using graphite and Kalrez^® seals—both seals experienced some problems during the experimental runs. The graphite seal leaked a very small amount and exhibited high friction on the axial piston, and the Kalrez o-ring seal decomposed to produce hydrofluoric acid, which caused some local vessel corrosion at the o-ring vessel wall contact.

Stress was applied to the rock sample by means of a 10-ton hydraulic jack (Enerpac C101) that moves the shaft against a titanium pedestal that supports the rock core. When the jack is activated through the use of a syringe pump (Teledyne Isco 260D), it applies force that compresses the rock core in the direction perpendicular to the asperity-propped fracture. The amount of sample displacement caused by fracture closing is monitored by a high-temperature linear voltage differential transducer (LVDT; Active Sensors LT-0951-010) that is attached to the shaft.

Fluid is introduced into the system from a reservoir of deionized water using a high-pressure syringe pump (Teledyne Isco 500D). A temperature equalizer bath is used to ensure that the fluid inlet and outlet temperatures are kept cool. Water is injected into the top hole of the core sample at a constant flowrate of 0.01 ml/min, except during system start up when the volume needs to be filled. Water exiting the margins of the rock core fracture moves out of the pressure vessel, through a cooling bath, through a back pressure regulator (BPR) that is controlled by a high-pressure syringe pump (Teledyne Isco 500D), and then collected for analysis using a fraction collector under ambient temperature and pressure (the reacted water is not recirculated through the system). Fluid inlet and outlet pressures are monitored continuously with pressure transducers, and low differential pressure (below ∼100 kPa) was independently measured.

The pressure vessel and its contents were heated by band heaters surrounding the pressure vessel, mediated by a ∼1.6 mm-thick copper shell to distribute the thermal flux. The entire vessel, including the heaters, was thermally insulated with 2.54 cm-thick glass wool having a thick aluminum foil surface, to stabilize the temperature of the system. A thermocouple and temperature controller were used to monitor and maintain a constant temperature for the experiments. This system was designed for temperatures as high as 300 °C; however, the actual experiments were conducted at 150 and 250 °C.

Four different rock samples were selected for these experiments—a granite, a metasedimentary mudstone, a rhyolite ash-flow tuff, and a silicified rhyolite tuff. The first sample was obtained from the Stripa mine in Sweden, the metasedimentary mudstone and the rhyolite ash-flow tuff were obtained from the BCH-03 well from the Bradys geothermal field, and the silicified rhyolite tuff was from the DP 35-13 well from the Desert Peak geothermal field. Table 1 provides additional details regarding the sample recovery depths and mineralogy based on petrographic study and more detailed mineralogic characterization from the literature for the same rock type. Figure 2 depicts photomicrographs for each of these samples.

Rock cores were prepared for each of the four starting rock samples. Before inducing a fracture, these cores had a diameter of 50.4–50.8 mm, and a length of 50.8–50.9 mm. First, a 2.6 mm diameter hole was drilled along the core axis up to about 60% of the length, to serve as an injection port for the fracture flow experiment. Subsequently, localized compressional force was applied at 12 points on the center circumference of the sample via ball bearings driven by set screws, to generate a tensile fracture perpendicular to the core axis. The core halves were then each jacketed with a passivated titanium ring (Fig. 3) with an O.D. of 63 mm and an I.D. closely matched to the rock cores, which was heated and then slid onto the core. The objective of these sleeves was to apply lateral constraints to the deformation of the core sample, such as to avoid development of new fractures in the vertical direction. The rim of this ring has several holes for a locking pin, which allow the core sections to be rotated and slightly misaligned at a designed angular offset during the experiment.

The two fractured core faces for each sample then were carefully scanned using an optical surface profilometer (Nanovea PS-50) to characterize the topography of the fracture surfaces. Because of the limitation of the area covered by a single scan, data were collected for overlapping quarters of the rock face at different scan heights, and processed using ImageJ software [30]; the data were then normalized and stacked using a matlab script, and the compiled quarter sections were stitched together to make a circular profile surface. The fracture surfaces from both top and bottom were registered to align them and subtracted from each other to estimate the fracture aperture; these scans could be used to estimate fracture apertures and asperities for rotated surfaces, and pre- and post-test fracture topographies could also be compared with each other to determine changes in the fracture surface resulting from the hydrothermal experiments.

The two halves of the titanium-jacketed core sample were put together with a slight rotational misalignment to create a self-propped fracture, then placed inside of the pressure vessel, with the drilled hole facing upward to facilitate fluid injection into the sample. To ensure the flow was going through the hole and not along the sample surface, a gold foil ring gasket was placed at the inlet, which was squeezed and created a chemically inert seal once the axial stress was applied. Prior to conducting the long-term hydrothermal experiments, an initial test under ambient conditions and increasing effective stress was conducted for each sample to identify the effects of initial fracture closure under compression and resulting changes in fracture permeability. Pore pressure (10.3 MPa) was applied to the samples under room temperature by injecting fluid into the sample chamber, and then the effective stress on the rock core was increased stepwise to a value of 20.7 MPa. This effective stress would be equivalent to a lithostatic load at 1408 m depth assuming a bulk rock density of 2.5 g/cm³ and a water density of 1.0 g/cm³, if the fracture was oriented horizontally. The change in mean fracture apertures of the samples with increasing effective normal stress was calculated by subtracting the measured closing displacement of the fractures (using the linear voltage differential transducer) from the initial fracture aperture of the uncompressed sample.

Due to sample compaction, decreases in the calculated mean fracture apertures were observed for all four samples (Fig. 4). The Stripa granite had the smallest initial aperture, because a smaller rotational mismatch than other samples was applied to this sample to create the aperture. Also note that the titanium sleeves applied to this sample lacked an alignment/locking pin that would have prevented the core halves from rotating back into alignment during the test. The rhyolite ash-flow tuff was loaded into the test cell twice—the sample was prematurely unloaded after a week due to a problem with the experimental apparatus, and the resulting initial aperture of this sample for the second test was about 35% lower than for when it was first tested.

Once the samples were fully placed under the uniaxial loading normal stress, the samples were heated and changes in inlet and outlet pressures and fluid chemistry were monitored over the duration of the experiments. Summaries of the experimental parameters for each of the tests that were run are listed in Table 2. All experiments were run with a pore pressure of 10.3 MPa, an effective normal stress of 20.7 MPa, and a fluid injection flowrate of 0.01 ml/min.

During the experiments, the vessel temperature, applied confining pressure, the sample displacement, and the differential pressure between the fluid inlet and outlet were measured continuously. Using the equations listed previously (Eqs. (2) and (3)), a fracture aperture (assuming a parallel plate model) and a wormhole tube radius were calculated based on the measured differential inlet and outlet pressures; these values as well as the measured properties mentioned above are plotted for each of the experimental runs in Fig. 5. These computed hydraulic apertures provide a convenient conceptual model. It is important to note that the computed hydraulic apertures from pressure differentials shown in Fig. 5 do not match the much larger average geometric apertures shown in Table 2. This is expected because the two are equivalent only for the case of flat plates with uniform separation, which is rarely observed in nature and was not the case for the experiments.

After the test runs were completed, the rock samples and flow-through fluids were characterized to evaluate the effects of water–rock interaction at elevated temperatures and mechanical compression. The topography of each of the fracture surfaces was rescanned with the profilometer to determine if any detectable changes caused by mineral solution and/or precipitation had occurred. The surfaces were also examined using scanning electron microscopy (SEM) to look for evidence of mineral dissolution and precipitation. The identification of mineral phases on these surfaces (Figs. 7–10) was achieved through the use of semi-quantitative energy-dispersive X-ray spectroscopy, crystal morphology, and X-ray diffraction (XRD) using powdered minerals collected from the surface of the fractures. These measurements were somewhat challenging due to the small amount of precipitated material, but when used in combination (XRD signal, chemical composition, and crystal morphology), allowed for reliable identification of mineral phases.

Aliquots of fluid effluent from each of the flow-through experiments were analyzed to determine changes in fluid chemistry caused by mineral dissolution and/or interaction with the experimental apparatus. An automated sample collection system was used to collect composite fluid chemistry samples every 6 hours—these samples were 3.6 ml in size. (In some cases, the sample collection time was extended to 12 h—the flowrate was still the same). The effluent samples were analyzed for cation contents using an ICP-MS after being acidified with 2% ultrapure nitric acid. Anion concentrations were analyzed using an ion chromatograph—samples were diluted 5–12 times using ultrapure deionized water prior to analysis. Deionized water used in the experiments as input water had no detectable anions, cations, or silica. Bicarbonate concentrations were not directly measured due to experimental limitations—instead, these were estimated by assuming the positive charge balance of the samples was compensated entirely by bicarbonate. A blank system sample was run at 250 °C for 88 h using a Ti plug instead of a rock core to estimate how the experimental apparatus may have affected the effluent fluid chemistry—the blank shows that elevated Al (31.9 mg/l—from the Grade 5 Ti pressure vessel), F (96.9 mg/l—from the elastomer seals), and SiO₂ (168 mg/l—likely from contamination from previous samples in the outlet line) did occur. Summaries for average water chemistry compositions are presented in Table 3, and representative plots of changes in effluent chemistry over time for each experiment are presented in Fig. 6.

There are some common trends in many of the effluent fluid composition plots depicted in Fig. 6. In most cases, there are initial rises in the concentrations of Na, K, Li, and Ca, followed by a decrease and a plateau that may represent a quasi-steady-state value. For the rhyolite ash-flow tuff sample (and also for the silicified rhyolite—not depicted), the elevated Cl and Na values observed in the effluent at the start of the run are interpreted to represent mobilized geothermal brine trapped as pore water—thus the Na and Cl compositions of these samples have been corrected to remove this to properly evaluate the dissolved mineral components (see Table 3). The short-lived spikes in K, Li, and Ca may in part also be related to the capture of residual geothermal brine from the rock pores or to the dissolution of highly soluble phases. The dissolved silica contents of the three samples heated to 250 °C are consistent with equilibrated silica contents using the quartz (no steam loss) geothermometer [32] for the experimental run temperature, which should yield a value of 467 ppm SiO₂. The Stripa granite sample, which was heated to 150 °C, has silica contents (60–70 mg/l) that are about ½ that corresponding to quartz solubility at that temperature (125 mg/l); this may reflect slower mineral-fluid equilibration rates at this lower temperature. The metasediment effluent exhibits unexpected compositional variability that was not observed in any of the other three experiments—the cause for this is not known.

Granite run results—As noted prior to starting the hydrothermal tests, the Stripa granite had the lowest initial fracture aperture, with a calculated hydraulic fracture aperture starting out less than one micron following loading of the sample: this value gradually decreased during the first day of the experiment (matched by an increase in the pressure differential between the inlet and outlet ports), and then appeared to stabilize over time. The Stripa granite sample had a lower run temperature (150 °C) and a shorter run duration (506 h), so the amount of hydrothermal reaction for this sample should be less extensive than the other three rock samples. The effects of lower temperature for this run are reflected in the much lower silica concentrations of the effluent (64.3 mg/l); the F concentrations are lower due to reduced thermal degradation of the elastomer seal, and lower corrected TDS values for the fluid (193 mg/l). Prevalent albite dissolution and diaspore precipitation are seen along the fracture surface (Fig. 7)—zeolites were not observed.

Rhyolite ash-flow tuff run results—This run had the longest duration of the experiments (3150 h); more than four times longer than any of the other runs. Significant dissolution of plagioclase (andesine/labradorite) is observed on the fracture surface from the dissolution images. Several phases were observed as precipitating minerals on the fracture surface: these include phillipsite (a calcium zeolite), illite/smectite, and Fe oxyhydroxides (Fig. 8). The elevated sodium and silica in the effluent water is consistent with plagioclase dissolution. While there was a systematic increase in the LVDT measurement (suggesting a compression of the sample by ∼100 µm), there was little change in the inlet and outlet pressure differential, and thus little variation in the computed fracture aperture during the experiment. Fluid flow may have occurred both within the fracture and in the secondary porosity generated by mineral dissolution (plagioclase, quartz, and calcite) within a ∼400-µm thick zone adjacent to the fracture surface.

Metasediment run results—The effluent from this experiment had the highest silica concentrations (482 mg/l), and also had the highest calculated bicarbonate concentration (247 mg/l), more than 60 mg/l higher than that of the next highest sample, which is consistent with dissolution features observed on calcite grains (Fig. 9). Precipitation of aluminosilicate phases such as phillipsite and chlorite also appears to have occurred on the fracture surface. There is a fairly extensive zone (∼150-micron thick) of mineral dissolution immediately adjacent to the fracture surface—this appears to have formed by selective dissolution of reactive mineral phases (calcite and quartz). Although the sample exhibited little variation in the LVDT measurement, apart from a systematic shift at ∼360 h, which would record a change in aperture caused by overall compression of the sample, there were repeated increases and decreases in differential pressure between the inlet and outlet ports, indicating fluctuating changes in fracture permeability during the run. These changes appear to have been caused by selective mineral dissolution occurring just below the fracture surface, which led to the development of dynamically evolving anastomosing flow paths.

Silicified rhyolite tuff run results—This sample experienced very little in the way of compaction over the ∼500 h duration of the experiment. A small but gradual difference between the inlet and outlet pressures was observed over the first half of the experiment; however, a sharp rise in pressure differential occurring during the last 250 h of the run reflects a marked decrease in fracture permeability, resulting in a decrease in computed fracture aperture from ∼200 µm during the early stages of the experiment down to a value of ∼25 µm at the end. SEM and BSE images (Fig. 10) show clear evidence of dissolution of plagioclase, calcite, and quartz, along with precipitation of the zeolites phillipsite and wairakite, together with illite/smectite and Fe oxy-hydroxides on the fracture surface.

Weathering and precipitation on the surface were not evenly distributed: preferential flow paths highlighted by the precipitation of new material were clearly present in the samples. SEM images were taken where the dissolution and precipitation events were more clearly visible.

A mass balance calculation using the dissolved species in the effluent can also be used to estimate the amount of net dissolution that has occurred. The amount of effluent produced during each of the runs (calculated by the flowrate of 0.6 ml/h and the experiment duration) along with the concentration of total dissolved solids can be used to estimate the total net amount of dissolution (total mass dissolution—mass precipitation). Table 4 lists the estimated total amount of dissolution for each of the experimental runs, using the corrected chemistry values in Table 3 that have removed the Al and F contributed by the experimental apparatus. The average value of equivalent fracture thickness that has been dissolved ranges from 11.4 µm (for the granite sample) up to 230 µm (for the rhyolite ash-flow tuff). This calculation is an upper estimate, as it assumes that all mineral dissolution occurs along the fracture surface, and does not occur along the 2.6 mm diameter hole drilled through the upper core section. It also assumes that mineral dissolution is distributed evenly along the fracture surface—however, it is most likely to be concentrated in the center of the fracture near the injection hole. The two samples with the highest calculated amounts of net dissolution (and the longest run times) also had the most secondary dissolution porosity observed in the backscatter electron images (Figs. 8 and 9). The ash-flow tuff had a ∼400 µm thick zone with fairly extensive dissolution—this presumably was present on both sides of the fracture. Using a total dissolution zone of 800 µm and the calculated dissolution thickness of 259 µm, this would then equate to a secondary dissolution porosity of ∼29%. The metasediment had a ∼150 µm thick zone with fairly extensive dissolution adjacent to each face of the fracture—using a total dissolution zone of 300 µ and the calculated dissolution thickness of 62 µm, this would then equate to a secondary dissolution porosity of ∼21%. The granite sample, which had the least amount of net dissolution based on the total dissolved solids in the effluent (11.4 µm average thickness), had very little evidence of dissolution beyond the fracture surface itself, but was also run at a significantly lower temperature (150 °C).

The samples run at 250 °C exhibited significantly more mineral dissolution than the granite sample, which was reacted at 150 °C. The main mineral phases that appear to have undergone dissolution in reaction with the deionized water were plagioclase, calcite, and quartz—the dissolution of these phases generated secondary porosity near and along the fracture surface. In the granite run, diaspore rosettes were observed—these may be in part precipitated due to the elevated Al concentrations of the fluid resulting from reaction with the Al-bearing titanium pressure vessel. For the higher temperature runs, the secondary minerals that were observed to have precipitated along the fracture surfaces include the Ca-zeolites phillipsite and wairakite, chlorite, illite/smectite, and Fe oxy-hydroxides.

There was no clear-cut correlation between sample compaction (as measured by the LVDT and the pressure differential between the inlet and outlet ports), which reflects changes in permeability for the fracture. This suggests that a pure parallel plate model, where the closing of the planar fracture aperture results in a reduction in fracture permeability, does not adequately reflect the observed changes. Instead, it appears that flow is more likely to be somewhat channelized along the fracture in areas with larger asperities, and also occurs in developing flow channels created by dissolution that results in a zone of enhanced secondary porosity immediately adjacent to the fracture surface. These zones were observed to be up to 400 µm in thickness on either side of the fracture. In general, fracture permeability decreased over time as a result of the combination of sample compaction resulting from the confining pressure, mineral dissolution, and mineral precipitation. There was fairly good correlation between increasing net mineral dissolution based on the amount of dissolved solids in the effluent and the observed thickness of the dissolution zones adjacent to the fracture surfaces within the samples—these amounts increased with increasing reaction temperatures and experimental run times.

The dissolved solids contents of the effluent fluids are considerably lower than those of equilibrated brines from the Bradys geothermal field (Table 3), from where two of the rock samples were obtained. This is due to the abbreviated time that the injected deionized water has to react with the host rock. There are two primary system volumes where fluid is in contact with the rock core: the injection hole drilled through the center of the core, and the fracture. The injection hole volume can be calculated from the hole length (∼25 mm) and diameter (2.6 mm)—this works out to ∼133 mm³. Using the injection rate of 0.01 ml/m, this volume would be filled in 13 min. The fracture volume is a bit more challenging to estimate, as the aperture varies significantly. Using an aperture value of 600 µm (fairly representative of most of the post-loaded samples) (Fig. 4), and the core diameters (50.4 mm), a fracture volume of 1200 mm³ is calculated, which would require 2 h to fill. Thus, the total time that a batch of injected water would be in contact with the rock would be on the order of 2.2 h—given that flowrates will vary within the fracture surface (especially for water within the secondary porosity adjacent to the fracture), these values will vary. Smaller aperture values would result in even shorter circulation times. Even with these short fluid residence times, the effluent appears to be close to equilibrium with respect to silica for sample runs conducted at 250 °C.

The conditions of the rhyolite ash-flow tuff experiment were used to constrain two different coupled process models that were simulated using the TREACTMECH simulator, which combines the parallel TOUGHREACT thermal-hydrological-chemical (THC) code [33,34] with geomechanics [35,36]. The code allows modeling tensile and shear failure within elements using tangential moduli in place of original moduli in elements of open failure. Tangential moduli are the derivatives of stress with respect to strain, after accounting for the effects of failure [37]. The open fracture between the core halves was considered as an initial open tensile fracture. The permeability was calculated locally at the mm to sub-mm scale accounting for the local aperture at this resolution. Therefore, the fracture-scale permeability includes a significant part of the roughness from sub-mm scale asperities. Furthermore, the large aperture changes from geomechanical and chemical effects have a much larger effect on permeability than roughness. The simulations were performed in two steps: thermal-hydrological-mechanical (THM) simulations to capture the geomechanical effects, and then THC simulations using the end results of the THM model as a starting point.

The first model simulated the THM effects of the experiment on the fracture. The reactive components of the experimental apparatus (including the core of rhyolite ash-flow tuff, the discrete fracture, flow-through elements, and titanium pressure vessel) were fully represented in the simulations. The resulting 36 × 36 × 12 element grid (Fig. 11) has 0.00508 m vertical spacing for most of the grid and 0.0108 m vertical spacing for the bottom three rows of elements. The grid blocks consist of three distinct elements: water flow zones, the titanium alloy piston and sleeves jacketing the core sample, and the rhyolite ash-flow tuff rock sample. The different properties assigned to each element are listed in Table 5.

The initial fracture aperture for the model was created by using the scanned unreacted fracture surfaces for the two core halves and then offsetting the matched cores (determined by best fit of the two halves) by a rotation of −2.3 deg. A three-point asperity contact model was used to calculate the initial fracture aperture. In the simulation, the sample was then subjected to a uniform fluid pressure of 10.3 MPa and a 10.3 MPa stress state (along with a small gravitational gradient) at 250 °C. Fluid flow was allowed from boundary cells around the outer edges of the fracture, where abovementioned pressure and temperature conditions were maintained, with no fluid flow boundary conditions for the other elements. The top boundary was held at a fixed position z, and a uniform normal stress was applied at the bottom of the system (corresponding to the titanium piston in the actual apparatus), ramping up from the initial 10.3 MPa average vertical stress to 31 MPa over a 20 s period (corresponding to an average vertical effective stress of 20.7 MPa), and then held at these conditions for another 20 s to allow conditions to stabilize. With increasing vertical stress, strain on the three asperity points propping open the fracture initiates fracture closure. When an element in the fracture is closed, the original unmodified elastic properties are retained, but with zero tensile strength. Initial geomechanical modeling using shear and tensile strengths based on measured unconstrained compressive strength for the rhyolite ash-flow tuff [25,37] resulted in shear failure of the three elements containing the contact point asperities that originally propped open the fracture. The model was modified to locally increase shear and tensile strengths of the tuff to prevent material failure beyond the initial tensile failure occurring in these fracture elements.

Resultant simulated fracture apertures under 20.7 MPa average effective stress reflect a significant decrease in fracture aperture, with 25 elements in contact along the fracture plane. At zero vertical effective stress, the effective aperture is 247 µm. With the application of 20.7 MPa vertical effective stress, the effective aperture decreases to 12 µm. Most of the model fracture surface elements with less than a 170 µm aperture decline are in contact, as illustrated in Fig. 12.

The results of the THM model were used as the starting point for the THC simulations, using the resulting 20.7 MPa aperture field depicted in Fig. 12. Thermodynamic and kinetic data were chosen from simulations of the Newberry enhanced geothermal system [39] and experiments on Desert Peak rhyolitic tuff. The deionized injection water was pre-equilibrated with CO₂ and O₂ at atmospheric pressure to capture the experimental conditions. The initial mineral assemblage consisted of calcite, dolomite, quartz, montmorillonite (Mg, Ca, Na, K), illite, clinochlore, plagioclase (albite-anorthite), sanidine, biotite (annite-phlogopite), apatite, and hematite—this was guided by XRD analyses reported by Ayling et al. [27] for a correlative sample from a depth of 1189 m from the same borehole. All minerals were allowed to dissolve/precipitate under kinetic constraints. The rhyolite core was given a porosity of 0.0345, and a permeability of 1.82 × 10⁻¹⁹ m². The initial pore water was assumed to have 0.1 M NaCl based on observations of elevated Cl⁻ concentrations in effluent over the early part of the experiment, consistent with typical reservoir salinities observed in geothermal reservoirs in this area (Table 3). The THC simulation was run for a period of one month (720 h)—the actual experiment was conducted over a longer period of time (3150 h). Figure 13 shows plots of mineral precipitation, dissolution, and calculated changes in fracture porosity and permeability. Fluid flux distributions show flow predominantly in the hole and in the main fracture, as expected. The two main phases observed to dissolve were plagioclase and quartz, with most dissolution occurring near the injection hole at the center of the fracture. Based on these simulations, laumontite (a zeolite very similar to phillipsite, but which has better thermodynamic data) also formed, consistent with the experimental observations. Large increases in porosity (>0.15) and permeability (up to two orders of magnitude) were observed near the injection hole at the center of the fracture, with much smaller changes seen as the fluid was dispersed radially across the fracture. These results indicate that the geochemical effects are concentrated in the center region of the fracture. Here, the fluid/rock ratio will be higher due to the radial flow and the fluid is least equilibrated with the rock by being closest to the inlet.

The simulated fluid chemistry for the rhyolite ash-flow tuff sample can be compared with the observed effluent compositions. For the simulations, the core pore water was assumed to have 0.1M NaCl, very similar to the reported geothermal brine composition for Desert Peak (Table 3). Figure 14 depicts the simulated and observed chemical trends for Na, K, SiO₂, and Cl. There is fairly good agreement for the long-term concentrations of Na and SiO₂ for these samples—the variations for the other components may be a result of not capturing all of the precipitating mineral phases in the simulations.

Several observations from the experimental and modeling work merit discussion, including comparing model system behavior versus experimental system behavior, changes in rock mechanical integrity, geochemistry of multi-mineralogic systems under gradients, changes in the fracture hydrological structure, and extension to the field.

Our geochemical modeling of the experimental system resulted in behavior similar to our experimental observations. Due to modeling constraints, it is not realistic to incorporate the full topology of the fracture surfaces at the measured resolution, nor the mineralogy at that scale, and the critical scale of the processes that are occurring may be different. In spite of these limitations, the model results and experiment results show important similarities. As expected, dissolution would be expected to be most prominent near the inlet, and indeed that was observed in both the model and experiment. The modeling indicated a decrease in plagioclase, quartz, and calcite near the inlet, and experimental observations showed the same. SEM images show increased porosity near the inlet, and the model shows the same thing. We were not able to experimentally quantify permeability at the model gridblock scale for comparison with the model, but it stands to reason that the more open porous structure resulting from mineral dissolution should increase in permeability over the initial case as is indicated in Fig. 13. The precipitation of laumonite was indicated by the numerical modeling, whereas phillipsite was indicated in the experiment. These two minerals are quite similar, lending confidence in the THC modeling. Even with the relatively short fluid residence times in the experimental apparatus (2.1 h) and the use of deionized water, most effluent samples from the higher temperature experiments appear to have achieved equilibrium with quartz.

Precipitation of secondary minerals still occurred even near the inlet in spite of the deionized water and low fluid residence time. In the field, injection of cooler water into fractures will result in disequilibrium between the mineral phases and the water. Even if equilibrated water could be introduced at the injection well for that temperature condition, the change in temperature away from the inlet will result in disequilibrium, and the condition at the injection well will change as the temperature changes. The complexity of the mineralogy and geochemistry will result in chemical gradients that could extend some distance into a fracture and will change over time as dissolution, precipitation, and temperature changes occur. Because these gradients and dissolution and precipitation kinetics can be captured in the modeling, the effects of such complex processes on fracture fluid flow can be evaluated [39,40].

Observations of increased porosity near the fracture surface due to dissolution will cause reductions in rock mechanical strength there. The more porous rock is expected to be far less rigid than the original solid rock. Imaging studies [41] have shown that host rock near asperities tends to be less geochemically impacted by dissolution and thus maintains its mechanical strength better. That is consistent with reduced flowrates in tighter apertures and asperity contacts as well. This could result in dissolution and some erosion at locations away from asperities, having the combined effect of increasing aperture and enhancing local flow, but possibly generating particles that might be filtered out downstream reducing flow there. The dissolution process tends to be self-enhancing to a point, as increased aperture and flow will lead to increased dissolution. Ultimately, the asperities will carry more stress, and will become weaker as dissolution occurs there, and they will break, resulting in a reduced aperture.

Precipitation of minerals and filtering of particles will tend to slow flow. Particle filtering will occur at narrow apertures, reducing the already lower flow there. Mineralization will limit or occlude flow wherever it occurs [16]. These processes affect the flow path as well. Water flowing through an EGS reservoir must be well-distributed to evenly extract heat from the rock. Self-enhancing dissolution could result in undesirable wormhole creation; this process would be limited by diffusive mass transfer from the location in the rock where dissolution is occurring to the flowing fracture, and relative flow velocities in the fracture, which are in turn affected by numerous processes. Precipitation and filtering will tend to alter flow paths as well; however, this will likely strengthen asperities. Because of these multiple processes, changes in the properties of interest for geothermal energy extraction of the resulting fracture are not immediately obvious and modeling is important to understand them on a larger scale. These models will need to take into consideration the specific water chemistry and rock mineralogy as well as reservoir pressure and temperature conditions.

Creating and sustaining fracture flow pathways is critical to the long-term performance of EGS reservoirs. Effective heat transfer from hot rock to a working fluid over long time periods requires creating (through stimulation) and retaining open fractures. While there have been numerous laboratory experiments and numerical modeling of field-scale processes focused on fracture stimulation for EGS [7,18–20,39,40,42–45], there has been relatively little focus on the long-term evolution of fracture permeability related to continued stress under hydrothermal conditions that could lead to changes in fracture aperture caused by compression, shear failure, mineral dissolution, and mineral precipitation. THC modeling of field tests conducted at the Soultz EGS site [46,47] suggest that mineral dissolution and precipitation can impact reservoir permeability and long-term flow behavior.

This study provides key insights by conducting hydrothermal-mechanical experiments on EGS candidate rock types at relevant pressures and temperatures to evaluate temporal changes in fracture permeability. All of the samples in our tests experienced net dissolution, with increased dissolution associated with longer run times and higher temperatures. The samples with the largest amount of dissolution exhibited patchy bands of dissolution where plagioclase, calcite, and quartz crystals were preferentially dissolved, forming a spongy region of secondary porosity between 150 and 400 µm thick along both sides of the fracture margin. The dissolution will reduce the mechanical strength of the fracture walls, allowing some closing of the aperture, and also result in the generation of fines. Some secondary mineral precipitation was also observed in the form of zeolites, clays, and Fe-oxyhydroxide phases. Such precipitation will also affect permeability by occluding flow at the formation locations or at downstream pinch points. THMC modeling revealed that most of the dissolution was focused near the inlet port at the center of the fracture, where the water/rock ratio was highest and the fluids were most out of equilibrium with the minerals. This dissolution led to localized increases in fracture porosity, but the experiments suggested that the fracture experienced an overall decrease in permeability over time resulting from the combined geomechanical and hydrothermal processes. The results of this combined experimental and modeling study can help inform the design and operation of field-scale EGS reservoirs to ensure that THMC processes do not result in premature degradation of the EGS fracture network.

Preliminary results of this project were presented at the Stanford Geothermal Workshop. We thank Jonathan Icenhower and Li Yang for their analytical work supporting this project. We thank Ormat Technologies, Inc. for providing the meta-sedimentary, rhyolite and silicified tuff core samples and associated rock composition information. This work was supported by the Assistant Secretary for Energy Efficiency and Renewable Energy (EERE), Office of Technology Development, Geothermal Technologies Program, of the U.S. Department of Energy under Contract No. DE-AC02-05CH11231.

There are no conflicts of interest

The datasets generated and supporting the findings of this article are obtainable from the corresponding author upon reasonable request. The data and information that support the findings of this article are freely available at.² The authors attest that all data for this study are included in the paper. Data provided by a third-party listed in Acknowledgment.