A wide, systematic variation of sedimentary geothermal gradients has been previously observed along the northern continental shelf of the Gulf of Mexico. From east to west, geothermal gradients change from 25 to 30 °C/km off Alabama to lower values (15–25 °C/km) off eastern Louisiana and to higher values (30–60 °C/km) off Texas. In order to assess the mechanism responsible for this variation, the present study first compiled an extensive bottom-hole temperature database from over 6000 wells in the northern continental shelf and constructed a more detailed geothermal gradient map than those published previously. Second, basin models were then constructed for three areas within the continental shelf (off Texas, Louisiana, and Alabama) that show differing geothermal gradients. A basin model is a mathematical model that simulates the heat transport through the crust and the sediments of a basin in the context of its geologic evolution. Previous researchers proposed two possible causes for the observed geothermal gradient variation in the northern continental shelf. The first was the thermal effect of sedimentation: areas with faster sediment accumulation result in low geothermal gradients, and vice versa. The second was that basal heat flow (heat flow that enters from the igneous crust to the bottom of the sediments) varied across the continental shelf. The present study finds that sedimentary geothermal gradients in these areas are primarily impacted by two competing mechanisms associated with sediment accumulation. One is the radiogenic heat production within the sediment that adds to the total heat budget upward through the sedimentary column. The other is the transient effect of fast sediment accumulation, which results in reduction in the upward heat flow. Off Louisiana, the transient effect prevails, and hence the area shows the lowest geothermal gradients. Off Texas, due to slower sedimentation, the positive contribution by radiogenic heat is most significant. Off Alabama, because the sediments there are not as thick, the overall contribution of radiogenic heat is less. The models show that the thermal effects of sedimentation are large enough to explain the observed variation in geothermal gradients. Therefore, corresponding variation in basal heat flow is not required.

The northern continental shelf of the Gulf of Mexico has been drilled extensively for petroleum exploration and production for the past several decades. One of the byproducts of the drilling activity is the tens of thousands of bottom-hole temperature (BHT) measurements obtained from these boreholes. The BHTs are the most abundant subsurface temperature measurements, and hence researchers have utilized them to infer the geothermal setting of this region (e.g., Bebout and Gutierrez, 1981; Blackwell and Richards, 2004; Forrest et al., 2005; Nagihara and Smith, 2008; Nagihara, 2010). However, these researchers have also recognized the limitation in the usefulness of BHTs. A BHT as measured is substantially below the pre-drilling formation temperature, because the circulation of drill fluid lowers the temperature around the wellbore while drilling. A BHT may easily be 10–20 °C lower than the pre-drilling formation temperature at the same position in a deep (>3 km) well (e.g., Harrison et al., 1983). Therefore, previous researchers have proposed various methodologies for estimating pre-drilling formation temperatures by applying some kind of correction to BHTs (see the review by Beardsmore and Cull, 2001). Another problem with use of BHTs is that until computer data entry became common in the mid-1990s, BHTs were often incorrectly reported in log headers.

Using only BHTs corrected for the effect of drill fluid circulation, Nagihara and Smith (2008) produced a map of geothermal gradients for the northern continental shelf of the Gulf of Mexico. These authors recognized two geographic trends. First, from east to west, geothermal gradient changes from 25 to 30 °C/km off Alabama to lower values of 15–25 °C/km off eastern Louisiana and to higher values of 30–60 °C/km off Texas (Fig. 1). Second, geothermal gradients tend to decrease toward the outer continental shelf. Based on these observations, these authors suggested that the regional variation trends in geothermal gradients are attributable primarily to the rapid and regionally variable sediment accumulation during the Cenozoic (e.g., Galloway et al., 2011). Sedimentary particles absorb geothermal heat and become equilibrated with the underlying particles, as they become buried deeper. Fast sediment accumulation causes retardation of this thermal equilibration process and results in a lower geothermal gradient through the pile of newly accumulated sediments than in that accumulated at a slower rate (e.g., Sharp and Dominico, 1976; Hutchison, 1985). Nagihara and Smith (2008) suggested that geothermal gradients are lowest off Louisiana mainly because, since the Miocene, sediment accumulation has been faster there than elsewhere in the northern continental shelf of the Gulf of Mexico. The geographic variation trends of the geothermal gradients match quite well with that of the thickness of Pleistocene sediments (Fig. 2). These authors, however, did not rule out a possibility that basal heat flow (i.e., the upward heat flow entering the sedimentary column from the igneous basement) and/or heat production associated with radioactive decay of uranium, thorium, and potassium within the sediments may also vary across the continental shelf and affect the sedimentary geothermal gradients.

Further understanding of the geothermal setting of the northern continental shelf of the Gulf of Mexico is important for a number of reasons, such as assessing risks of thermogenic hydrogen sulfide gas (Nagihara and Smith, 2005), sandstone reservoir quality degradation due to quartz cementation (e.g., Taylor et al., 2010), and potential for geothermal power generation (Griggs, 2005; John et al., 2006). The present study has two goals: (1) to characterize the regional trends in geothermal gradients of the northern continental shelf in more detail, using the BHT database recently compiled for the National Geothermal Data System (Nagihara et al., 2013), which is substantially larger than the one used by Nagihara and Smith (2008); and (2) to examine whether or not the thermal effects of sediment accumulation alone can explain the magnitude of the observed regional variation in geothermal gradient by constructing vertical heat transport models at three locations of different geothermal gradients in the continental shelf.

The National Geothermal Data System (NGDS) is a worldwide-web portal of geothermal data, recently established by the United States (U.S.) Department of Energy (geothermaldata.org). The present authors led the effort to compile BHTs from the U.S. Exclusive Economic Zone of the Gulf of Mexico in 2011 through 2013 and delivered them to the NGDS node managed by Southern Methodist University (Nagihara et al., 2013). For the Gulf of Mexico data set, BHTs were tabulated from the log headers of over 6000 wells distributed in the northern continental shelf and slope of the Gulf of Mexico (Fig. 1). The logs were obtained from the U.S. Bureau of Safety and Environmental Enforcement (BSEE).

The so-called “Horner plot” method was used for our estimation of pre-drilling formation temperatures. This method is widely used by geothermal researchers, and it has already been described in detail elsewhere (Lachenbruch and Brewer, 1959; Deming and Chapman, 1988; Beardsmore and Cull, 2001). Several other BHT correction methods have been previously proposed, and a few of them are considered more reliable than the Horner plot method (Hermanrud et al., 1990; Beardsmore and Cull, 2001). However, these methods require data that are not usually reported in well log headers such as the thermal properties of wellbore sediments, and therefore their use is limited (see review by Nagihara et al., 2014).

The Horner plot method requires that BHTs be measured multiple times at one depth in a well. Such a series of measurements should show BHT gradually increasing with time and recovering toward its pre-drilling state. Mathematical extrapolation of the BHT recovery trend to infinite time yields an estimate of the pre-drilling formation temperature, also known as the virgin rock temperature (VRT). Not every BHT can be corrected in this manner, because the vast majority of the wells examined here did not report BHTs multiple times at one depth. For the present study, we used VRT estimates from ∼1590 wells located within the northern continental shelf (Fig. 1). The volume of the VRT database used in the present study is more than four times that of Nagihara and Smith (2008).

In obtaining local geothermal gradients, we used the same methodology used by Nagihara and Smith (2008). Because of the requirement of multiple BHTs measured at one depth, one well typically yielded only one or two VRT estimates at different depths. Such a small sample size was not adequate for reliable determination of the geothermal gradient at that location. Therefore, we combined VRTs from several wells located within close proximity to one another into a group (typically 5- to 20-km radius). Such a group of wells together yielded ∼10 VRTs, and that allowed determination of the geothermal gradient of that vicinity. In the map of Fig. 2, we show three examples of such groups of wells that collectively define the geothermal gradients for their vicinities (Fig. 3). This methodology worked relatively well in the continental shelf, where seafloor is flat and sedimentary thermal properties do not vary drastically within a short horizontal distance. Grouping of the wells was carried out in such a way that the VRTs within each group yielded an easily recognizable temperature versus depth trend. In some of these well groups, spacing of the wells can be ∼50-km radius or more, because of the sparsity of the well data. For such groups, the well-grouping approach may introduce smoothing effects in the final map.

Some quality control measures were applied to the data. First, for some of the wells, BHT values did not change at all for tens of hours after they were shut. We viewed such cases as erroneous reporting of BHTs and excluded them from the present study. It is widely accepted that, for deep (>1000-m) wells, drilling fluid circulation cools the formation around the wellbore (e.g., Jaeger, 1961; Wooley, 1980; Shope et al., 2012). Therefore, after the well has been shut-in, BHT should begin to increase as the wellbore formation gradually recovers to its original thermal state. Second, grouping the BHTs from wells in close proximity to one another helped us identify obvious outliers in their temperature versus depth trends. Such outliers may have resulted from incorrectly reporting of BHT measurements or the limitation of the correction method (Nagihara et al., 2014).

The local geothermal gradients were spatially interpolated to yield a continuous map for the northern continental shelf as shown in Figures 1 and 2. The natural-neighbor interpolation (Kennedy, 2009) was performed on the ArcGIS software by the Environmental Systems Research Institute (ESRI). The two regional trends initially observed by Nagihara and Smith (2008) can be seen. From east to west, geothermal gradient values are intermediate (25–30 °C/km) off Alabama-Mississippi, low (15–25 °C/km) off eastern Louisiana, and high (30–60 °C/km) off Texas. Off Louisiana and east Texas, where the continental shelf is widest, geothermal gradients decrease toward the outer shelf. Smaller-scale variations that were not obvious in the 2008 map can also be seen. Even though geothermal gradients are generally high in the Texas shelf, the highest values tend to occur along the Corsair growth fault zone (Fig. 2), as previously suggested by Nagihara (2010). The new map shows that this trend continues farther east to the shelf off western Louisiana.

Farther east, there are several bull’s-eyes of high geothermal gradient off Louisiana. Given that there are a number of salt structures in this area, it is possible that these high-temperature anomalies are associated with geothermal heat flow funneled through high-conductivity salt bodies (e.g., O’Brien and Lerche, 1984; Nagihara, 2003). However, it is difficult to determine whether or not these anomalies exactly coincide with particular salt structures, based on the previously published salt maps alone (e.g., Ewing, 1991; Diegel et al., 1995; Watkins et al., 1996).

Origin of the high geothermal gradients along the Corsair fault zone has been previously discussed but not fully understood. The Wilcox growth fault zone on the Texas coastal plain, which runs parallel to the Corsair growth fault zone, shows a similar high-gradient trend (Blackwell and Richards, 2004). There, it is very likely that hot fluids, expelled from deep, over-pressured sediments, migrate upward along some of the faults and locally elevate the sedimentary temperature (McKenna, 1997). A similar mechanism is possible for the Corsair fault zone (Taylor and Land, 1996). Alternatively, Rangin et al. (2008) propose a failed rifting event that took place concurrently with the initiation of the Corsair growth-faulting in the Oligocene to early Miocene, and they suggest that the high geothermal gradients are the remnant from the rifting that may have affected the sediments in the entire Texas continental shelf. Such a hypothesis may be tested by thermal modeling through the crust and the sediments of the area.

Nagihara and Smith (2008) suggested that the east-west variation trend in geothermal gradient across the northern continental shelf is primarily due to the corresponding variation in sediment accumulation rate in the Cenozoic. The faster the accumulation rate is the lower the geothermal gradient. In sedimentary basins such as the Gulf of Mexico, three processes primarily control sedimentary temperature distribution on such a regional scale. The first is the thermal transient effect of rapidly accumulating sediments. The second is the basal heat flow, which is influenced primarily by the tectonic process that created the basin. The third is the radiogenic heat production within the sediment.

In order to further understand the origin of the observed geographic variation in geothermal gradients, we have constructed computational models of the heat transport through the sediments that incorporate these three processes (so-called “basin models”) for three areas in the northern continental shelf (Fig. 2) representing the high geothermal gradients off Texas (Area 1), the low geothermal gradients off Louisiana (Area 2), and the medium-range geothermal gradient off Alabama (Area 3). A basin model simulates the heat transport up the sedimentary layers as they accumulate over time. Building a reliable basin model requires knowledge of the tectonic process that created the basin, the sediment accumulation history from the beginning, the thermo-physical properties of the sediments, and the present-day temperature distribution within (or heat flow through) the sediment.

The Gulf of Mexico basin formed during the breakup of the supercontinent of Pangea. Rifting began in the Late Triassic. Following the deposition of massive evaporites in the Middle Jurassic, a brief period of seafloor spreading took place in the Late Jurassic (e.g., Buffler, 1991; Salvador, 1991; Hudec et al., 2013). The structural framework of the igneous basement of the basin has changed little since then. The present-day northern continental shelf lies over the continental crust stretched during the initial rifting (Marton and Buffler, 1994). During the Cretaceous, extensive carbonate platforms were constructed along the margins of the basin. In the Cenozoic, sedimentary influx from the Mississippi drainage basin has dominated sedimentation in the northern Gulf margin. Depocenters gradually shifted from west to east along the Texas-Louisiana shorelines (Woodbury et al., 1973; Galloway et al., 2011). Since the Miocene, the continental shelf off Louisiana has seen fastest accumulation of clastic sediments. There, the total thickness of the sediment reaches 15 km (Peel et al., 1995; Radovich et al., 2011).

In the Texas-Louisiana continental shelf, where no boreholes are deep enough to reach the igneous basement, availability of seismic-reflection sections imaging deep sediments that allow interpretation of their ages is crucial for constraining the basin models. We chose Areas 1 and 2 (Fig. 2) mainly because of the availability of previously published seismic sections that show the oldest sediments (McDonnell et al., 2009 for Area 1 and Radovich et al., 2011 for Area 2). In Area 3, a number of boreholes reach the Jurassic sections (e.g., the Norphlet and the Smackover formations), which allow us to rely less on seismic information.

We performed the basin modeling on the software package PetroMod 1D of Schlumberger. The software uses the finite element method in obtaining numerical solutions to the one-dimensional heat conduction equation with a source term (e.g., Hantschel and Kauerauf, 2009). Temperature distribution through the sediment column is obtained at each time step as sediment accumulates at the top. Using the exponential compaction curves (e.g., Sclater and Christie, 1980), the model accounts for sediment compaction and resultant change in its thermal and physical properties. The software package comes with a library of rock thermal properties and their dependence to temperature and porosity (Hantschel and Kauerauf, 2009). The sedimentary rocks in the study area are divided into two general groups based on the findings from previous studies (e.g., Galloway et al., 2011): the Mesozoic sections are rich in carbonates, while the Cenozoic sections are dominated by clastic sediments. For each sedimentary layer in the model, we chose a parameter set most appropriate for the expected lithology and observed porosity.

The PetroMod 1D software is mainly designed to quantify the transient thermal effects of sediment accumulation and the radiogenic heat production within sediments (e.g., Hutchison, 1985). It does not account for the upward advective heat flow associated with dewatering of the compacting sediment. Previous studies suggest that the compaction-driven fluid flow does not contribute significantly to the total budget of the heat flow through the sedimentary column in the continental margin setting such as the present study area (e.g., Harrison and Summa, 1991).

The rate of sediment accumulation, its change through time, and sediment compaction significantly influence the heat flow through the sedimentary column (e.g., Hutchison, 1985; Hantschel and Kauerauf, 2009). For the present study, the thicknesses and the ages of the sediments were determined from the wireline logs from the wells in these areas and the fossil records from the core samples compiled by the BSEE (www.bsee.gov). Porosity and lithology of the sediments were also obtained from wireline logs. For deeper sediments not reached by boreholes, their stratigraphy and ages were obtained from the previously interpreted seismic sections (Peel et al., 1995; Radovich et al., 2011). The software estimates thermal conductivities and heat production rates of the sediments from the porosity and the lithology of the sediments, using the previously proposed models correlating them (Hantschel and Kauerauf, 2009). The estimated properties are built into the rock properties library of the software.

The PetroMod 1D software uses the “pure-shear” lithospheric stretching model (McKenzie, 1978) to obtain the basal heat flow. In our study areas, rifting took place in the Late Triassic through the Late Jurassic (ca. 220 Ma to ca. 160 Ma). The spike in basal heat flow associated with the rifting has almost completely subsided since that time. Radiogenic heat production within the stretched, igneous crustal layer is now the main controlling factor of the basal heat flow. We refer to Sawyer et al. (1991) for the present-day crustal thickness.

We constrain the basin models primarily by matching the overall trend of the model-predicted porosity-versus-depth relation against that of the log-derived porosities and matching the model-predicted present-day temperatures against the VRTs estimated from the NGDS database (Fig. 4). It should be noted that the aim of these models is to compare the three areas in the continental shelf that have large differences in their overall geothermal setting (Fig. 2). These models are not intended to reproduce small-scale variations such as a porosity kick associated with a particular over-pressured formation at a certain depth. In addition, because the aforementioned model parameters have some uncertainties, the models presented here should be viewed as good representation of the overall geothermal setting of our study areas, not a precise, unique representation of the reality.

In assessing how sediment accumulation, basal heat flow, and radiogenic heat production affect present-day sedimentary temperature distribution, it is helpful to track the heat-flow budget upward through the entire sedimentary column. For example, the model for Area 1 (Fig. 5) shows that the basal heat flow at the bottom of the sedimentary column is ∼44 mW/m2. Heat flow is near-constant through the Mesozoic sections. Through the Tertiary sections, heat flow gradually increases, because radiogenic heat produced within these sediments is added to the total heat budget. In the northern continental shelf of the Gulf, the Cenozoic sediments are dominated by clastic sediments, while the Mesozoic sediments are more dominated by carbonates (e.g., Galloway et al., 2011). Shale produces more radiogenic heat than most other sedimentary rocks (e.g., Beardsmore and Cull, 2001). Therefore, on the average, clastic sedimentary sections are more heat-producing than carbonate sediments. Farther upward in the upper Miocene through the Pleistocene sections, heat flow no longer increases, because (1) these sediments are more porous and thus produce less radiogenic heat per volume than the deeper, more compacted sediments; and (2) a greater portion of the heat budget is used for thermally equilibrating the more rapidly accumulating sedimentary particles.

In summary, as geothermal heat enters from the basement and travels upward through the sedimentary column, some heat is added by radiogenic heat production, and some heat is lost by the transient effect of sedimentation. If sediment accumulation is slow, the positive contribution by the radiogenic heat production prevails and results in a net gain of heat flow upward. Fast sediment accumulation results in a net loss. As shown in Figure 5, the model for Area 1 shows a net gain of nearly 9 mW/m2. In contrast, the model for Area 2 shows that all of the heat gained by radiogenic production is cancelled by the rapid sediment accumulation in the Miocene through the present. The result is a net loss of ∼14 mW/m2 in the total heat-flow budget. Area 3, where sediment accumulation has been much slower throughout the Cenozoic, sees a net gain of only 4 mW/m2.

Figure 5 also shows that the three areas have similar basal heat-flow values of 42–44 mW/m2. This suggests that the east-west variation in geothermal gradients observed in the northern continental shelf can be accounted for almost entirely by the balance of the two competing thermal effects of sedimentation. Geothermal gradient is lowest off Louisiana because the heat loss by sedimentation is much greater than the heat gained by radiogenic heat production. Between the two areas off Texas and Alabama, the former gains more radiogenic heat. As seen in Figure 5, the total thickness of the Cenozoic (primary clastic) sections in the area off Texas is nearly three times as thick as that in the area off Alabama. Even though the Texas shelf experienced faster sediment accumulation in the Quaternary and lost some of the heat gained at greater depth, it yields greater net gain relative to the Alabama shelf.

Nagihara and Smith (2008) first reported the geothermal gradient variation along the northern continental shelf of the Gulf of Mexico with 25–30 °C/km off Alabama to lower values of 15–25 °C/km off eastern Louisiana to higher values of 30–60 °C/km off Texas. To better understand possible causes of the observed variation, we first compiled an extensive database of bottom-hole temperatures from 6000+ wells in the northern continental shelf. We applied corrections for the disturbance associated with drilling to BHTs from ∼1600 wells and updated the geothermal gradient map of Nagihara and Smith (2008). Second, we built basin models for three areas in the northern continental shelf, off Texas, Louisiana, and Alabama. These models simulate the thermal evolution of their sediments and reveal that the two competing mechanisms that affect the upward heat flow through the sedimentary layers are present in all of the three areas. The first is that the fast sediment sedimentation in the Cenozoic reduced the upward heat flow through the sedimentary column. The second is that the thick sedimentary layers contribute a significant amount of radiogenic heat to the total heat budget. For the area off Louisiana, even though it has the thickest sediment pile, the negative effect of the fast accumulation completely cancels out the radiogenic heat added. For the area off Texas, sediment accumulation in the Cenozoic has been considerably slower, and thus the total heat budget shows a net gain of heat; therefore, geothermal gradients are highest there. For the area off Alabama, even though sediment accumulation has been slowest of the three in the Cenozoic, because the heat-producing sediment pile there is not as thick as that in the area off Texas, radiogenic heat production does not add as much heat. Therefore, the observed east-west variation in geothermal gradients across the continental shelf does not require corresponding variation in basal heat flow.

This work was supported by a grant from the U.S. Department of Energy. Comments by two anonymous reviewers and Dr. John Holbrook, the associate editor, were very helpful in the revision of this manuscript. Discussion with Dr. Will Gosnold was helpful in the interpretation of the basin model results.