Seismic properties in conductive and convective hot and super-hot geothermal systems

Seismic methods contribute to the exploration of geothermal areas and characterization of existing geothermal resources. Seismic velocity and attenuation depend on the pressure and temperature conditions of the geothermal systems, which are closely related to the properties of the rock frame and geothermal fluids. We calculate the seismic velocities and attenuation in terms of the subsurface distribution of the confining and pore pressures and temperature, assuming that the heat transfer from below is convective or conductive. The pore pressure is assumed hydrostatic. In hydrothermal systems the temperature is calculated assuming the boiling point condition at the specific pore pressure down to the reservoir. Beneath the reservoir it is assumed constant in convectively heated systems and following a constant temperature gradient in conductively heated systems. In Enhanced Geothermal Systems (EGS) conductive heat transfer and constant temperature gradient are assumed. We present three application examples, considering simplified subsurface models to describe the geothermal systems beneath the production wells. The seismic wave properties are calculated using the rock's mechanical Burgers model and the Arrhenius equation to take into account rock-properties-variability with temperature and the Gassmann model for fluid saturating the porous rocks.


Introduction
Seismic methods are used for the exploration of geothermal areas and characterization of existing natural and potentially enhanced geothermal resources (EGS) (Majer et al., 2007;Aqui and Zarrouk, 2011;Iovenitti et al., 2013) basically, to understand and verify the conceptual model, explore new resources, and identify and characterize the fracture/fault systems (e.g., Majer et al., 2007;Hashida et al., 2001). Appraisal methods make use of surface and borehole seismic measurements (Poletto et al., 2011;Niitsuma et al., 1999) for the location and characterization of geothermal reservoirs. Passive seismic methods are also useful, e.g., in volcanic areas (Blanck et al., 2016), even though they are less accurate for well siting. Simiyu (2010) applied microseismic methods to geothermal exploration, Bannister et al. (2015) obtained the 3D seismic velocity and attenuation to investigate the deep geothermal resources in a volcanic area. Surface reflection seismic response together with passive seismic has been used for the characterization of deep structures and shallow layers in different geothermal regions (e.g., Batini et al., 1983;Bannister, 1992;Majer, 1978;Majer and McEvilly, 1982;Henrys et al., 1986;Kristinsdóttir et al., 2010).
Knowledge of the physical properties and microstructure of the geothermal rocks is important to understand their seismic response as a function of pressure, temperature and fluid properties . To support deep geothermal exploration, extensive works have been conducted on numerical simulations of geothermal systems focused on developing improved geophysical models to describe the properties associated to the presence of brittle-ductile transition (BDT) (e.g., Montesi, 2007) or supercritical fluids Reinsch et al., 2017;Farina et al., 2016Farina et al., , 2017Hashida et al., 2001). Carcione and Poletto (2013) modeled the BDT using the Burgers mechanical model to describe the ductility effects, and the octahedral stress criterion and the Arrhenius equation to calculate the flow viscosity as function of temperature and pressure. Carcione et al. (2014) developed an algorithm to simulate the full-waveform propagation including the effects of ductility and the temperature dependence of flow viscosity. More recently, Carcione et al. (2017) extended the theory and the simulation algorithm to poro-viscoelastic media using the Gassmann equation to take into account the presence of geothermal fluids.
On the basis of these studies, Poletto et al. (2018) calculated the sensitivity of seismic wave propagation to physical properties and https://doi.org/10.1016/j.geothermics.2019.05.005 Received 4 February 2019; Accepted 7 May 2019 temperature variation using the Burgers model augmented with the Arrhenius relation and integrated in a modified Gassmann model to determine the sensitivity of the elastic properties, stiffnesses, impedance and attenuation to temperature, including frequency dependent effects related to permeability (Manning and Ingebritsen, 1999), fluid mobility (Batzle et al., 2006;Zhubayev et al., 2013), and squirt flow (Carcione et al., 2018a,b;Carcione and Gurevich, 2011;Gurevich et al., 2010). The analysis of the seismic quantities is extended to the dependence of elastic moduli to temperature variations (Jaya et al., 2010;Poletto et al., 2018). All these aspects are related, with different relevance depending on the specific geological context, to the seismic characterization of a geothermal reservoir. From the conceptual-system point of view, a key issue in the characterization of a geothermal reservoir is to understand the nature, whether conductive or convective, of the heat-transfer mechanism. This nature, under different conditions in different geological scenarios, determines the fluid-rock temperature regimes at depth and may have important impacts for the exploitation of geothermal reservoirs and for the production of geothermal resources (Edwards et al., 1982). In this work, we use the numerical approach of Poletto et al. (2018) to calculate seismic velocity and attenuation in poro-viscoelastic media where the temperature profiles are calculated considering models of conductive and convective heat flow systems. The aim of this analysis is to clarify a series of key issues, namely: 1. The seismic characterization and identification of the conductive and convective processes. 2. Determination of the most sensitive visco-elastic properties. 3. The influence of super critical fluids in very hot geothermal regions. 4. How the seismic properties are affected by melting.
The potential use of the proposed approach is studied by using known petrophysical rock properties, and envisaged for different geothermal sites and temperature-pressure profiles, including existing and potential Mexican geothermal sites investigated in the framework of the GEMex project (GEMex, 2016).
We first review the poro-visco elastic seismic modeling theory, extended to include temperature, pressure as well as melting effects. Then, we introduce the characteristic properties of conductive and convective geothermal systems, provide examples of seismic analysis with known crustal rock parameters, and discuss the results for different reservoir scenarios assuming pure water as geothermal fluid.

Theory: the Burgers-Gassmann model
The anelastic behaviour due to shear deformation and plastic flow is described using the Burgers mechanical model (Carcione and Poletto, 2013;Carcione, 2014;Carcione et al., 2014), whose frequency dependent shear modulus is where ω is the angular frequency, = − i ( 1), μ 0 is the relaxed shear modulus of the Zener element describing the brittle material, τ σ and τ ϵ are seismic relaxation times expressed by the relations where Q 0 is the minimum quality factor and τ 0 is a relaxation time such that ω 0 = 1/τ 0 is the center frequency of the relaxation peak. The shear viscosity η S , which describes the medium ductile behaviour, is obtained by the Arrhenius equation (e.g., Carcione et al., 2014;Montesi, 2007) accounting for thermodynamic rehological effects. This quantity is related to the steady-state creep rate ϵ˙by where σ o is the octahedral stress (e.g., Gangi, 1981Gangi, , 1983Carcione et al., 2006;Carcione and Poletto, 2013), A ∞ and n are constants, E is the activation energy, R G = 8.3144 J/mol/K is the gas constant and T is the absolute temperature. All these parameters govern the ductile creep and melting behaviour. The octahedral stress is where the σ's are the stress components in the principal system. In our analysis we assume that these components correspond to the vertical (v) confining stress, and the maximum (H) and minimum (h) horizontal tectonic stresses, given by where x and y are the horizontal coordinates, z is depth, z = 0 corresponding to the surface, ρ is the medium density, g is the acceleration of gravity, and, omitting for simplicity the dependence on lateral coordinates x and y, The parameter ξ = ξ(x, y, z) ≤ 1 accounts for the additional effects due to tectonic stresses (Carcione and Poletto, 2013) and ν = ν(x, y, z) is the Poisson's ratio of the formation. In addition to the Burgers mechanical model, the Gassmann model allows us to predict the low-frequency limit of the wet-rock bulk modulus with complete saturation . This model considers the elastic effects due to the mineral components, rock-frame bulk modulus (K s ), dry-rock elastic moduli (K m and μ m ), porosity (ϕ) and pore fluid bulk modulus (K f ). According to the Gassmann model, the saturated-rock bulk and shear moduli are given by (e.g., Carcione, 2014) The elastic behaviour of granular material composing the dry rock depends on pressure and is non-linear. In general, velocity rises with increasing confining pressure, and this effect is probably due to pores/ cracks closure and hence increases the density. At low effective pressure, cracks are open and easily closed with an increase in pressure (Chen et al., 2015). The elastic moduli typically vary as power functions of mean stress (Houlsby et al., 2005). The pressure dependence of dryrock-moduli can be expressed as (Kaselow and Shapiro, 2004) where K 0 and μ 0 are the dry-rock bulk and shear moduli at infinite (i.e., very high) confining pressure, a 1 , a 2 , p 1 and p 2 are constants and p d = p c − p o is the differential pressure, p o and p c are the pore and the confining pressure, respectively. In Eq. (11) we use μ 0 = μ B , where μ B = μ B (ω) is the frequency dependent Burgers shear modulus given by Eq.
(1) that takes into account the rock modulus variation in the presence of high temperature and melting. The complex seismic velocities are given by is the bulk density, and ρ s and ρ f are the grain and fluid densities, respectively. The phase velocities and quality factors are respectively. These equations have been used, together with terms including frequency-dependent effects related to permeability and fluid flow, to study the sensitivity of seismic properties to temperature variations in a geothermal reservoir (Poletto et al., 2018), assuming linear gradient temperature models typical of conductive environments. In this work, we extend the analysis to other realistic geothermal models.

Conductive and convective systems
Temperature and presence of fluids are key conditions for the characterization of a geothermal reservoir. In particular, temperature is a key parameter for the evaluation of geothermal resources. The temperature profile versus depth can be determined only with direct access to the rocks, through measurements in drilled wells. In the absence of direct measurements, it is possible to use temperature models derived from geophysical, geological and geochemical prospecting (e.g., Manzella, 2010). The models can be characterized by different behaviours in different environments. An important classification is based on the distinction between conductive and convective heat-transport systems summarized below (Axelsson and Steingrímsson, 2012;Moeck and Beardsmore, 2014).

Conductive systems
In a purely conductive system, the heat flow remains almost constant with depth, as stated by the first law of thermodynamics, while the thermal gradient varies according to the conductive properties of the rocks (Beardsmore and Cooper, 2009). This means that the temperature profiles varies with depth. EGS are cases of conductive systems, as in situ permeability is too small to allow the movement of fluids.

Convective systems
A convective system differs from a conductive one, because in its upflow zone the fluid specific enthalpy is nearly constant. For a liquidonly rising fluid this implies a nearly constant temperature with depth. The same is approximately valid for a vapor-only rising fluid if the pressure within the reservoir is nearly constant. For a two-phase rising fluid (liquid and vapor), where boiling occurs as the pressure decreases upwards, the temperature and pore pressure are interrelated following the boiling point with depth curve (BPD) down to the geothermal reservoir, which is the exploitation target. This is due to the convectivefluid exchange mechanism. Convection cells are emanated by a deeper heat source, which usually is a cooling magma chamber located beneath the geothermal reservoir. The heat transfer mechanism from below can be either convective, where convection cells extend down to the brittle-ductile transition (BDT) level (convectively heated hydrothermal system), or conductive (conductively heated hydrothermal system).
Convective geothermal systems can have different behaviours depending on the fluid regimes, phase and pressure. We distinguish the following subcases (Donaldson, 1982): Liquid-dominated system.The deeper reservoir fluid is in the liquid phase with temperature below the critical point of water, even though the wells may deliver a two-phase fluid. An indicative graph of fluid properties shown on the Mollier pressure enthalpy diagram (Elders and Fridleifsson, 2010), as it ascends from the deep heat source to the surface is presented in Fig. 1, assuming isenthalpic upward flow. Vapor-dominated system. The deeper reservoir fluid is in the vapor phase, while the wells deliver mainly steam vapor. Two indicative graphs of fluid properties shown on a pressure enthalpy diagram, as it ascends from the deep heat source to the surface are presented in Fig. 1, distinguished by the heat transfer mechanisms (convective or conductive).
The efforts in geothermal exploration of natural resources are typically focused on up-flow zones of hydrothermal convection system. The reference lowest bound of the temperature profiles is the linear thermal conduction curve with an average continental geothermal gradient of 30°C/km (Suzuki et al., 2014). The highest temperature profiles in those zones are normally limited by the BPD temperature curve (Haas, 1971). To simulate the upper part of the hydrothermal system, we calculate the temperature T(z) versus depth z in the conductive and convective models from surface to 2.25 km depth using the BPD curve (Henley et al., 1984), by the James equation (James, 1970), assumed as valid in the depth range 0.030 < z < 3.0 km,°= Below 2.25 km, we assume a temperature gradient of T = 120°C/km to reach the condition of T = 400°C at 2.7 km. In the deeper part, the temperature in the convective model is constant, whereas in the conductive model it increases with depth with the above temperature gradient. Below 2.25 km, the fluid properties as a function of depth are calculated every Δz = 100 m intervals from top downwards as follows. The fluid density is calculated from the temperature at depth z and the pressure of the overlying interval, while the pressure at depth z is calculated from the fluid density assuming fluid-static conditions, i.e., for The rock pressure is calculated from the rock density assuming lithostatic conditions, while the rock temperature is assumed equal to the fluid temperature. Fluid convective systems require fluid circulation conditions, related, to some extent, to permeability (Sorey, 1978;Cathles et al., 1997;Lipsey, 2014;Lipsey et al., 2016). Conversely, finite permeability is not necessarily required for a conductive system. Complex micro fracture systems, faulting and melting have a big influence (e.g., Saemundsson, 2013;Arnórsson, 2014). These factors may affect seismic wave propagation, such as variations in the Poisson ratio, attenuation and anisotropy associated to fracture orientation, that, in principle, can be investigated with seismic methods (e.g., Carcione, 2014). The effects related to temperature, including porosity and permeability, are discussed in Carcione et al. (2018a,b). Here, we model porosity, but we do not include permeability in the calculation of the seismic properties, even if the assumption of a convective system requires to consider permeability. We mainly focus on the effects induced by the thermodynamic properties, related to temperature and pressure changes expected for convective or conductive systems.
In the following analysis, we denote 'conductively heated hydrothermal systems' as 'conductive', and 'convectively heated hydrothermal systems' as 'convective'. We also consider the EGS scenario. Fig. 2 shows (a) the pore pressure and (b) the temperature profiles as a function of depth for the conductive (red lines) and convective (blue lines) models, respectively. In the first part, until depth of 2.25 km, the heat transfer is dominated by the convective mode. In the interval 2.25-2.7 km, we assume a transition from the BPD curve to 400°C by conduction, using the temperature gradient of 120°C/km. In the deeper part, we assume conductive and convective conditions. Both, these hydrothermal systems are 'vapor dominated'. The conductively heated hydrothermal system is characterized by much higher temperatures at depth, and by lower pressures, resulting from the lower fluid density. Depending on pressure and fluid properties, supercritical conditions may develop at these high temperatures . This implies different physical conditions that affect the seismic properties.

Examples
To characterize the convective and conductive systems seismically, we need the temperature conditions, the pressure of the saturating fluids, and the seismic properties of the formations, in terms of the elastic moduli of the dry rock and their dependence on the effective pressure. All these factors are involved in the Arrhenius equation (Eq. (3)). The characterization requires the collection of integrated information from geophysical, geological, and laboratory data.
Here, we focus on three representative examples, two of which related to a superhot geothermal reservoir and one related to an EGS geothermal reservoir. The properties of the geothermal fluids are assumed those of pure water. For a superhot geothermal system (examples 1 and 2), we consider the pore pressure and temperature conditions shown in Fig. 2, where below 2.7 km depth, we assume conductive and convective models.
We estimate the fluid properties by using the pore pressure and temperature profiles shown in Fig. 2, that characterize a vapor dominated superhot geothermal reservoir. We derive the density, the acoustic velocity and the bulk modulus of the fluid using CoolProp codes (Bell et al., 2014) based on the thermo-physical database provided by the National Institute of Standards and Technology (NIST). Fig. 3 shows the curves of (a) density, (b) acoustic velocity and (c) bulk modulus of the pore fluid for the conductive (red line) and convective (blue dashed line) models. The fluid properties are different only in the deeper layer, where the pressure and temperature conditions of the geothermal reservoir differ. Fig. 3d shows the difference between the fluid bulk moduli calculated for the two heat-transport modes (CV) , where subscripts 'CD' and 'CV' denote conductive and convective, respectively. For an EGS system (example 3) we consider a conductive-only model.

Example 1: crustal rock geophysical and thermal parameters
We consider a schematic model of a superhot geothermal reservoir with three layers. We use the geophysical and rheological properties of three crustal rocks reported in the literature including the dependence of the elastic dry-rock moduli on the differential pressure calculated in laboratory (Brace, 1965;Simmons and Brace, 1965;Popp and Kern, 1994). The Arrhenius parameters are derived from generic crustal formations (Fernández and Ranalli, 1997). The aim is to evaluate if, how and when the conceptual conductive and convective models can be seismically characterized and identified. We vary the porosity of the deeper layer to study the influence of the geothermal pore fluid in the presence of convective and conductive heat-transport mechanisms.
Then, we calculate the seismic properties of the saturated formation, assuming a given porosity, using the temperature and pore pressure conditions of a convective-liquid-dominated geothermal reservoir and compare the results with those of the conductive model. We vary also the Arrhenius parameters of the deeper layer, to see how the seismic properties change when the thermal properties melt the rock.
The profiles as functions of depth of the density, bulk modulus of the solid and porosity of the formations are shown in Figs. 4a, b and c, respectively. The dry-rock bulk and shear moduli at infinite confining pressure, the solid density and bulk modulus, with the sample's intrinsic porosity in brackets, are reported in Table 1. Properties for samples S-1 and S-2 (granite from Georgia and granite from Rhode Island) are reported in Brace (1965) and Simmons and Brace (1965), and those of samples KTB 61C9b, in Popp and Kern (1994). The frame bulk modulus at infinite confining pressure is calculated with the equation proposed by Krief et al. (1990), which relates the grain property to the dry-rock bulk modulus and the porosity considering the intrinsic porosity, i.e., neglecting compliant porosity (e.g., Poletto et al., 2018), of each sample and the dry-rock bulk modulus at infinite confining pressure. Constants a 1 , a 2 , p 1 and p 2 , used in Eqs. (10) and (11) to calculate the dry-rock bulk and shear moduli dependence on differential pressure, are reported in Table 2 for each sample. These constants are obtained from Brace (1965) and Simmons and Brace (1965) for samples S-1 and S-2, and from Popp and Kern (1994) for sample KTB 61C9b. The parameters that appear in the Arrhenius equation, used to characterize the three layers, are reported in Table 3, retrieved from representative rheological values for the wet upper crust reported in Fernández and Ranalli (1997). Castro et al. (2008) propose the shear seismic loss parameters for the crust in Southern Italy, where f is the frequency, an equation valid till 10 Hz. We calculate the relaxation times using Q 0 = 122 as the minimum quality factor, corresponding to a frequency f = 3 Hz. We use the Gassmann equation (7) to obtain the bulk moduli of the saturated formation (Fig. 5a). These moduli are substituted in Eqs. (15) and (16) to calculate the seismic properties as a function of depth, pressure and temperature. The difference of the saturated-rock moduli (CV) , is shown in Fig. 5b.
The compressional velocities of the conductive (red line) and convective (blue dashed line) geothermal systems are shown in Fig. 6a.
(CV) (Fig. 6b), has a maximum of 15 m/s, when the temperature difference between the two models is the highest. Similar trends can be observed for the compressional elastic moduli of the two heat flow models (Fig. 6c) and their difference (Fig. 6d). The shear velocities are shown in Fig. 7a, where the difference has a maximum of 9 m/s at the highest temperature gap (see Fig. 7b). The behaviour of the shear elastic moduli are shown in Fig. 7c and d.
The compressional (Q P ) and shear (Q S ) quality factors as a function of depth for both models are shown in Fig. 8a and b. The variability of the seismic properties is solely due to the properties of the fluid, since there is no melting. The fluid is in a vapor phase for both models, and the fluid properties do not significantly change in the deeper part, even if the temperature difference reaches 400°C. This results in small variations of the seismic velocities and quality factors. To investigate the variations of the visco-elastic quantities in a superhot geothermal system due to the presence of geothermal fluids in a vapor phase, we focus only on the deeper part below 2.7 km. Here, the pore pressure and temperature conditions change according to the conductive and convective heat transport mechanisms. This zone is modeled with the properties of the rock sample KTB 61C9b (Table 1). We vary the average porosity of the medium from 5% to 50%. Fig. 9a-c shows the difference in the bulk density (Δρ = ρ (CD) − ρ (CV) ), the difference in the compressional phase velocity ( (CV) ) and that of the shear phase velocity ( (CV) ). As the porosity increases, the effect due to the presence of the saturating fluid is more relevant, as can be seen in the variations of the bulk density and seismic velocities. The seismic quality factors are hardly affected.
For comparison, we calculate the seismic properties of a geothermal reservoir considering the temperature and pressure of a convective liquid-dominated (LD) system, where the temperature increases following the boiling-point to depth (BPD) curve until 1 km, where it reaches 300°C and then it remains constant. The petrophysical properties do not change, while the geothermal fluid properties and the dry-rock moduli variation with differential pressure are calculated using the pore pressure and temperature profiles shown in Fig. 10a and b, respectively. The fluid density, acoustic velocity and bulk modulus for the conductive (red line) and convective LD (blue dashed line) systems are shown in Fig. 11a, b and c, respectively. The difference between the fluid bulk moduli calculated for the two heat transport models ΔK f = K f(CD) − K f(LD) , where 'CD' and 'LD' denote conductive and convective liquid-dominated, respectively, is shown in Fig. 11d. This difference is approximately three orders of magnitude larger than that calculated for the vapor dominated reservoir (Fig. 3d). The change from vapor to liquid, related to the two different pressure-temperature conditions of Fig. 10, causes significant variations in the seismic properties. The compressional and shear velocities calculated with the conductive and convective LD mechanisms and their difference are    Table 3 Values of the Arrhenius parameters for the model layers.
Layer A ∞ (MPa) −n s −1 n E (kJ/mol) 1 1 0 −2 1.8 151 2 2 × 1 0 −4 1.9 134 3 2.9 × 10 −3 1.8 150 shown in Fig. 12. The maximum difference between the velocities of the conductive system, where the fluid is in a vapor phase, and the convective liquid-dominated system, is 246 m/s (P wave) and 137 m/s (S wave). The corresponding compressional and shear quality factors are shown in Fig. 13. The difference is small, with the compressional quality factor more sensitive.
To investigate the changes due to variations of the rock properties, we analyze the deeper layer, considering the pressure and temperature conditions of the vapor-dominated convective and conductive models. We vary only the Arrhenius parameters of the layer. We use the four sets of thermodynamic parameters reported in Table 4. Set A1 is the same used to calculate the seismic properties shown in Figs. 6-8 below 2.25 km depth. Set A2 is obtained from Violay et al. (2012). As for set A1, it characterizes rocks that melt at very high temperature, higher than 800°C. Rocks with the Arrhenius parameters and activation energy of set A3 (Poletto et al., 2018) and A4 (Fernández and Ranalli, 1997) start melting at around 700°C and 500°C, respectively. Fig. 14 shows the viscosity variation as a function of depth for the conductive (bold line) and convective (dashed line) heat flow models. The viscosity decreases with increasing temperature. Viscosity values lower than 10 11 Pa s are obtained for temperatures higher than 400°C in the conductive case, when the rock is characterized by sets A3 and A4. These low viscosity values are associated to the presence of melted material (Mavko, 1980;Solomon, 1972;Poletto et al., 2018), and this significantly affects the velocity (Fig. 15) and attenuation (Fig. 16) profiles. Conversely, in the convective model there is no melting because the maximum temperature is 400°C, and the rock viscosity does not change, remaining higher than 10 11 Pa s, without appreciable variations in the elastic moduli of the saturated rocks. Therefore, the seismic velocities calculated for saturated rocks with different thermal properties do not change.
More in detail, the compressional phase velocities calculated with the conductive and convective models are shown in Fig. 15a and b, respectively and their difference ΔV P in Fig. 15c. For rocks with Arrhenius parameters A1 and A2, there is no melt and the compressional seismic velocity increases by about 15 m/s when temperature increases from 400°C to 800°C. This variation is interpreted as due to the variations in the fluid properties. On the other hand, for rocks with properties A3 and A4, the compressional velocity decreases with increasing temperature with a maximum reduction of 1230 m/s and 1780 m/s, respectively, as expected in the presence of partially molten rocks (Solomon, 1972;Williams and Garnero, 1996;Carcione and Poletto, 2013;Poletto et al., 2018). The shear phase velocities calculated with the conductive and convective models are shown in Fig. 15d and e, respectively and their difference ΔV S in Fig. 15f. The difference is about 8 m/s for rocks with Arrhenius parameters A1 and A2. A large reduction in the shear velocity, reaching zero with complete melting, is observable in rocks with thermal parameter A3 and A4, when the conductive model is assumed, because the temperature exceeds the values at which melting starts. In principle, assuming the conductive  heat transport mechanism, the seismic observations can reveal the temperature variations and the presence of molten phases related to the sharp discontinuities in velocity and quality factor (Spetzler and Anderson, 1968;Carcione and Poletto, 2013;Poletto et al., 2018). Fig. 16 shows the calculated compressional and shear quality factors in the conductive and convective models and their difference. P-wave attenuation results from the relaxation of both the shear and the bulk moduli . The compressional quality factors (Q P ) calculated with the four sets are shown in Fig. 16a and b for the conductive and convective models, respectively, and their difference is shown in Fig. 16c. The shear wave attenuation results from relaxation of the shear modulus. Fig. 16d and e shows the shear quality factors (Q S ) of the conductive and convective systems, respectively, and their differences are shown in Fig. 16f. Seismic attenuations of rocks which start melting at temperatures greater than 800°C (A1 and A2) do not show relevant variations between the conductive and convective mechanisms. The attenuation, which is more related to thermally activated relaxation processes (Solomon, 1972) shows a very rapid decrease as the melting point is approached (Spetzler and Anderson, 1968;Carcione and Poletto, 2013), before reaching the brittle-ductile transition (BDT) zone, and this behaviour can be observed also in the convective model for set A4, which starts melting at about 500°C ( Fig. 16b  and d). Fig. 17 shows the seismic velocities and quality factors as function of temperature for the conductive model. The compressional and shear velocity curves ( Fig. 17a and b, respectively) abruptly decrease around 700°C and 500°C for rocks with parameters A3 and A4, respectively. The seismic compressional and shear attenuations start increasing at 400°C and 350°C for the two sets.

Example 2: Los Humeros geothermal system
We focus the analysis on the superhot geothermal field of the Los Humeros volcanic complex (Fig. 18) (GEMex, 2016), which is the largest active caldera located in the northernmost part of the eastern sector of the Trans-Mexican volcanic belt (Carrasco-Núñez et al., 2017). We consider the simplified lithological model of Gutiérrez-Negrin and Izquierdo-Montalvo (2010), derive the solid and dry-rock properties from the literature in addition to the average porosity and density (Aragón-Aguilar et al., 2017;García-Estrada, 1992).
Los Humeros is one of the oldest producing geothermal fields in Mexico (Arzate et al., 2018), the second after Los Azufres to produce electricity in the area of the Mexican Volcanic Belt (Prol-Ledesma, 1998). Many studies have been done to understand better the behaviour of the geothermal reservoir focusing on geophysics and geology (e.g., Arzate et al., 2018;Carrasco-Núñez et al., 2017;Urban and Lermo, 2013;Gutiérrez-Negrin and Izquierdo-Montalvo, 2010;Lermo et al., 2008;Cedillo Rodriguez, 2000), hydrogeology and hydrodynamic (e.g., Cedillo Rodriguez, 2000;Portugal et al., 2012), petrology and volcanology (e.g., Contreras et al., 1990;Ferriz and Mahood, 1984;Carrasco-Núñez et al., 2012), thermal and pressure conditions (e.g., Arellano et al., 2008Arellano et al., , 2000Verma et al., 2011;Verma, 1985). Arellano et al. (2000) studied the distribution of pressure and temperature of the Los Humeros geothermal field analyzing information from 42 wells drilled in the field. They proposed the existence of at least two reservoirs. The shallower one, located at 1.6-1.025 km above sea level (a.s.l.) is liquid-dominant with a pressure profile corresponding to a 300-330°C boiling water column. Arellano et al. (2003) used well pressure logs and observed that the data show high correlation with the boiling point pressure for this depth. The deeper one,    located at 0.85-0.1 km a.s.l., with low-liquid-saturation, has a temperature ranging between 300 and 400°C (Urban and Lermo, 2013).
Gutiérrez-Negrin and Izquierdo-Montalvo (2010) analyzed several Los Humeros wells, and showed that their temperatures follow the boiling point to depth curve down to a depth of about 2.5 km. However, beneath a depth of about 1.75 km, a cluster of production wells shows temperature values with a gradient of approximately 120°C/km, higher than those of the BPD curve and reaching a maximum temperature of 400°C at about 2.25 km depth. Pulido (2008) reported a maximum bottom-hole temperature of 395.4°C for well H-43, García-Gutiérrez   Table 4 Arrhenius thermodynamic parameter-sets used to evaluate the melting effect.
In this example, we assume the conductive and convective heat-flow models. For the geological model of the Los Humeros geothermal field, we consider the four main lithological units proposed by Gutiérrez-Negrin and Izquierdo-Montalvo (2010), according with the rock cuttings provided by the geothermal wells and based on previous works (e.g., Viggiano and Robles, 1988). These units are shown in Table 5. Starting from this lithological partition, we built a simple 1D four-layers model. For the calculation of the shear velocity, we use the reference ratio V P /V S = 1.76 (Lermo et al., 2008). We consider that additional petrophysical data, provided by the ongoing GEMex project (GEMex, 2016), will be used for a refinement of this investigation. Average porosity values have been assigned using porosity measurements in core samples (Aragón-Aguilar et al., 2017;Contreras et al., 1990).
The 1D velocities, density and porosity profiles are shown in Fig. 19 and the averaged seismic and thermodynamic properties are reported in Table 6. To select the corresponding Arrhenius parameters, we used the characteristic values of crustal formations, provided by Fernández and Ranalli (1997) for the first three layers. To evaluate the seismic response near the BDT, assuming the proximity of a magma chamber as a possible scenario, we consider the two sets of parameters ALH1 and ALH2 for the last and fourth layers, which correspond to two different behaviours of the rock at high temperatures. The set ALH1 characterizes a rock that melts at temperatures greater than 900°C (Violay et al., 2012), and set ALH2 characterizes a rock that melts at temperature around 700°C (Carcione and Poletto, 2013;Poletto et al., 2018).
The P-velocity values are taken from the literature values. Starting from these values, we calculate the frame bulk modulus K s , and using the porosity reported in Table 6, we invert the Krief relation given by Eq. (19) (Krief et al., 1990) to derive the bulk modulus of the dry matrix at infinite pressure (K 0 ). Using the shear velocity and the density, we calculate the shear modulus at infinite pressure μ 0 . The values of the tectonic parameter ξ = 0.8 in Eq. (6), of the minimum quality factor Q 0 = 122 and the frequency f = 10 Hz, used to calculate the relaxation times (Eq. (2)), are the same without variations in the different lithologies of the model. With these parameters, we obtain the seismic velocities and attenuations.
The seismic compressional and shear phase velocities of the saturated media corresponding to the conductive (red line) and convective (blue dashed line) case, for the model with the deeper layer characterized by set ALH1, are shown in Fig. 20a and b, respectively. The compressional and shear quality-factors are shown in Fig. 20c and d,  Table 4) is common to all panels.  Table 4) is common to all panels.  Table 4) is common to all panels.
respectively. In this case the variability in the seismic properties due to the different pressure and temperature profiles is only related to the fluid properties variations. At the maximum temperature and pressure differences, it is about 0.2% for the seismic phase velocities and 2% for  (Arzate et al., 2018). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)   Table 6 Seismic and thermodynamic average properties of the rock frame composing the media in the Los Humeros layered model. In the last layer we distinguish two thermodynamic rock behaviours related to two different Arrhenius parameters (A ∞ , n and E) sets named A1 and A2.  Farina, et al. Geothermics 82 (2019) 16-33 the seismic quality factors. The seismic velocities for set ALH2 are shown in Fig. 21a and b, respectively, while the compressional and shear quality factors are shown in Fig. 21c and d, respectively. In this case, the variability is mainly due to melting, i.e., close to the BDT zone. The phase velocities start decreasing at about 5 km depth, where the thermodynamic conditions allow rock melting. The maximum seismic phase velocities difference between the two analyzed temperature-pressure conditions is 20%. As expected, the quality factors related to the thermally activated relaxation processes (Poletto et al., 2018) show a rapid decrease starting approximately at 3 km, before reaching the melting point.

Example 3: Acoculco geothermal system
We consider the Acoculco Mexican site, where the geothermal conditions for a potential EGS site are investigated (Fig. 22) (GEMex, 2016). This area, situated near the town of Chignahuapan in the Mexican state of Puebla, is located in a volcanic complex that extends over the Ouebla-Hidalgo state boundary, in the eastern part of the Mexican Volcanic Belt (Canet et al., 2015;Pulido et al., 2010). Studies on two wells drilled in the Acoculco area have shown that the temperature profiles are linear, indicative of a conductive thermal regime (López-Hernández et al., 2009). We calculate the seismic properties of a fourlayer stratified model, assuming a conductive heat-transport mechanism.
The geothermal area of Acoculco is hosted by a volcanic caldera complex in the eastern part of the Mexican Volcanic Belt. Studies on this geothermal area have been performed to assess the feasibility of developing it as an enhanced geothermal system (EGS) for power generation (Pulido et al., 2010(Pulido et al., , 2011GEMex, 2016). The heat source is interpreted as related to the presence of magma, which heats the surrounding formation (Pulido et al., 2010). The Comisión Federal de Electricidad (CFE) drilled two exploratory wells in the southernmost area, and located not too far from each other, well EAC-1 in 1995 and well EAC-2 in 2008 reaching a depth of 1810 and 1900 m, respectively (Canet et al., 2015;Viggiano-Guerra et al., 2011). The area is characterized by the presence of active gas emissions. Thermal logs from the exploratory wells EAC-1 and EAC-2 show a conductive heat transfer  regime with bottom-hole temperature greater than 300°C (López-Hernández et al., 2009;Viggiano-Guerra et al., 2011).
In this example, we use the conductive heat-transport mechanism to calculate the temperature and pressure conditions shown in Fig. 23a and b, respectively. The temperature profile shown in Fig. 23a is that of well EAC-1 (blue bullets), obtained from Pulido et al. (2010). This profile is extrapolated in depth with a temperature gradient of 156°C/ km (orange line), which is the average gradient required to reach the bottom-hole temperature. The pressure in Fig. 23b is that of well EAC-1 (blue bullets, obtained from Viggiano-Guerra et al. (2011)). This pressure profile is in agreement with the pressure corresponding to a hydrostatic column of water (Viggiano-Guerra et al., 2011), and it is extrapolated in depth (orange line) calculating the hydrostatic pressure from water densities associated to the temperature profile of Fig. 23a derived from the NIST database. The water density, velocity and bulk modulus are shown in Fig. 24a, b and c, respectively. We can see that at about 2.4 km depth, the conditions are those of a transition from liquid (blue line) to supercritical (green line) water phase.
For the geological model, we consider a simplified 1D four-layers model representing the main lithologic units penetrated by well EAC-1 (López-Hernández et al., 2009) shown in Table 7. We assume compressional velocity and density values obtained from the literature. For the calculation of the shear velocity, we use the same ratio V P / V S = 1.76 used for Los Humeros. The 1D rock-frame compressional (blue line) and shear (red line) velocities and density profiles are shown in Fig. 25a and b, respectively. In the future, we will use data provided by the ongoing GEMex project (GEMex, 2016). We assume an average porosity of 6%, as proposed by Pan et al. (2016) for the Acoculco geothermal area. For the Arrhenius parameters, we use the same values of Los Humeros. The geophysical and thermal parameters are    Table 8. For the last layer, we consider two thermodynamic sets, AAC1 and AAC2, which correspond to two different behaviour of the rock at high temperatures. The Arrhenius parameters AAC1 = ALH1 characterize a rock which melts at temperatures higher than 900°C (Violay et al., 2012), and those of set AAC2 = ALH2 characterize a rock which melts at temperature around 700°C (Carcione and Poletto, 2013;Poletto et al., 2018). Fig. 26 shows the seismic compressional (a) and shear (b) phase velocities, and the compressional (c) and shear (d) quality factors, when the thermodynamic parameters do not allow melting (set AAC1, blue line), and when they allow the last layer to melt (set AAC2, red line). The seismic velocities of the medium, which melts at about 700°C, start decreasing at about 4 km depth, where we could expect the presence of the BDT. The effects on the quality factors start at a shallower depth. Recently, Calcagno et al. (2018) estimated the thermal gradient in the Acoculco area and the depth of the BDT zone at about 4 km depth below ground level. In this case, seismic measurements could in principle confirm this estimation.

Discussion
The paper presents a comprehensive analysis of the temperature effects on the seismic properties of geothermal reservoirs with different heat-transport mechanisms, focusing on convective and conductive systems in hot and very-hot regions. The analysis uses pure water as geothermal fluid, and the thermodynamic properties predicted by the Arrhenius equation. The study is based on simulation equations developed in previous works. One of the difficulties is to determine the input parameters for the simulations, related to the geological, geophysical, thermodynamic and fluid-transport properties. To afford this task, we have assumed a simplified layer system, using typical crustal parameters and basing on known scenarios of two reservoirs along the Mexican Volcanic Belt, one of them in operation.
We simplify the heat-recharge system, but the method can be useful to study more complex recharge reservoirs mechanisms, which are subjected to controversial investigations. Here, we assume equilibrium between the thermal properties of the rock frame and of the saturating fluid. Laboratory experiments are required to characterize the rocks near the BDT, where partial melt occurs. In this sense, this work is intended as a first step for future and more extensive characterizations of high temperature (HT) geothermal reservoirs, such as super-hot and EGS ones.

Conclusions
Our study presents a seismic characterization of convective and conductive geothermal reservoirs, with different thermodynamic properties dictated by the Arrhenius equation. The aim is to evaluate the influence of the geothermal mechanisms and temperature on the seismic properties, namely, seismic velocities, stiffness moduli and quality factors. The objective is also to discriminate between the two reservoirs at least in the hotter part, below the boiling point. The differences in the seismic properties are small when there is no melting, and are due to variations of the fluid properties. However, remarkable differences can be observed when passing from a vapor-dominated system to a liquid dominated system. Melt significantly affect the properties of the conductive reservoir, since in this case the temperature increases linearly with depth and highly affects the shear rigidity of the rocks. Conversely, in convective reservoirs, the temperature is constant with depth in the deepest region, and only partial melting can be observed for certain thermodynamic conditions.

Acknowledgments
This project has received funding from the European Union's Horizon 2020 research and innovation programme under grant agreement No. 727550.