Urey ratio and the structure and evolution of Earth's mantle

2.1.1. Bulk Silicate Earth and Its Derivatives

[9] The central concept in any type of global mass balance argument is the bulk silicate Earth (BSE), which is also known as the primitive mantle (PM). BSE refers to the solid Earth excluding the core, and it may be regarded as the primordial mantle formed after core formation but before the extraction of continental crust. By definition, all of present-day silicate reservoirs have been derived from BSE, so they must add up to BSE. We can readily recognize two major silicate reservoirs: the continental crust (CC) and the depleted mantle (DM), which have long been considered to be complementary to each other [e.g., Hofmann, 1988, 2003]. That is, the continental crust is the product of mantle melting, whereas the depleted mantle is the solid residue. A number of trace elements are incompatible with solid phases, so they tend to be highly concentrated in a melt phase, which is the continental crust in the context here. The mass of continental crust is only ∼0.6% of that of the mantle, but the continental crust is highly enriched in incompatible elements, so its contribution to mass balance is substantial. The depleted mantle is depleted in those incompatible elements with respect to BSE but not depleted in compatible elements such as major elements.

[10] Note that although the oceanic crust (OC) is also generated by mantle melting, it is usually not considered as an independent chemical reservoir in terms of global mass balance because it is quickly recycled back to the mantle by subduction. At a normal condition, suboceanic mantle beneath mid-ocean ridges has a potential temperature (hypothetical temperature of mantle adiabatically brought to the surface without melting) of ∼1350°C and starts to melt at the depth of ∼80 km, resulting in an ∼6-km-thick oceanic crust and an ∼70-km-thick depleted mantle lithosphere (DML) [e.g., McKenzie and Bickle, 1988; Langmuir et al., 1992]. DM in global mass balance refers to the source mantle before this melting takes place, and DML is even more depleted than DM. By mass conservation, the sum of OC and DML should be equivalent to DM if the effects of degassing and seawater alteration are taken into account.

2.1.2. Missing Heat Source Paradox and Its Solutions

[11] Among a few global mass balance arguments, the missing heat source paradox is most relevant to the thermal budget, so it is explained in some detail here. The continental crust is estimated to have 1.3 ppm U, 5.6 ppm Th, and 1.5% K [Rudnick and Gao, 2003], the heat production of which is ∼7.5 TW (assuming the continental mass of 2.3 × 10²² kg). According to the model of Salters and Stracke [2004], the abundance of those elements in the depleted mantle is 4.7 ppb U, 13.7 ppb Th, and 60 ppm K, and the corresponding heat production is ∼4.2 TW. Combined, these two reservoirs yield ∼12 TW. However, the bulk silicate Earth is commonly assumed to have 20 ppb U, 80 ppb Th, and 240 ppm K [e.g., McDonough and Sun, 1995], which give rise to the heat generation of ∼20 TW. Assuming that these values are reasonably accurate, the global mass balance does not seem to hold for these heat-producing elements. This is called the missing heat source paradox, and a key question is the following: where did we make a mistake?

[12] This failure of mass balance has been used to argue against whole mantle convection and to support the presence of a hidden reservoir in the deep mantle. This line of logic is understandable because the composition of the depleted mantle is based largely on the chemistry of mid-ocean ridge basalts (MORB). As mentioned in 3, mid-ocean ridge magmatism is the result of mantle melting at shallow depths (<100 km), so a composition model for the depleted mantle is fundamentally biased to the upper portion of the mantle. The above mass balance thus implicitly assumes that the mantle is well mixed by whole mantle convection so the shallow mantle sampling is not an issue. As we just saw, however, this assumption does not lead to a self-consistent mass balance, so it seems that we should abandon the notion of well-mixed mantle and call for a deep reservoir that is more enriched in the heat-producing elements. As the classical layered-mantle model gradually lost its popularity, various alternatives have been proposed for this “hidden” reservoir. The reservoir may occupy a substantial portion of the lower mantle, and its interface with the upper layer is somehow seismologically invisible [e.g., Kellogg et al., 1999]. Or the reservoir may be distributed as a number of blobs [e.g., Manga, 1996; Becker et al., 1999; Helffrich and Wood, 2001], though it is uncertain why such blobs do not contribute much to mid-ocean ridge magmatism (to an extent that enriched MORB becomes much more common than observed). Or it may correspond to the thin, D″ layer at the base of the mantle [e.g., Coltice and Ricard, 2002; Tolstikhin and Hofmann, 2005], but the layer might be too enriched in heat-producing elements to have been gravitationally stable in the past.

[13] The very first step we should have taken instead may be to examine the accuracy of the mass balance calculation that caused the paradox. As discussed in the following, every component in this mass balance has substantial uncertainty. First of all, the most recent compilation for continental crust composition by Rudnick and Gao [2003] assigns 30% uncertainty for trace elements, so the heat production in the continental crust is in the range of 5 to 10 TW. If one looks more closely at how the average composition of continental crust is estimated, even greater uncertainty may be more realistic. The continental crust is typically divided into three layers such as upper, middle, and lower crust, and all layers exhibit considerable regional variability. Global averaging is thus a difficult task. The 30% uncertainty adopted by Rudnick and Gao [2003] is for the upper crustal composition, which may be estimated most reliably because a large number of data are available. The deeper layers are expected to have larger uncertainties, which are left unquantified because of the paucity of data. One may argue that such uncertainty may not concern us because the deeper layers are often considered to be depleted in heat-producing elements. The conventional wisdom of heat production exponentially decreasing with depths [Lachenbruch, 1970], however, is just a model and does not seem to be supported by observations [e.g., Jaupart and Mareschal, 2003]. Also, intracrustal differentiation, which is likely to have occurred in the evolution of continental crust [e.g., Rudnick, 1995], may deplete the lower crust of K but probably not of U and Th, because these two elements are highly concentrated in accessory minerals such as zircon, which tend to be unaffected by crustal melting [e.g., Watson and Harrison, 1984]. Existing estimates for the lower crustal composition vary by an order of magnitude in terms of trace elements, and the “best estimate” according to Rudnick and Gao [2003] is based on crustal xenolith data, though they also acknowledge that it is unclear how well xenolith compositions could represent average lower crust and also that xenolith data are spatially limited. Given the potential of the continental lower crust therefore, the total heat production of continental crust may even be higher than 10 TW. Mantle geotherms based on xenolith data [e.g., Finnerty and Boyd, 1987; Russell and Kopylova, 1999] may be used to bound crustal heat production, but such constraints are not available globally. Though the heat production should not exceed the observed continental heat flow, this upper bound itself has been revised from 12 TW [Pollack et al., 1993] to 14 TW [Jaupart et al., 2007] on the basis of more recent heat flow data. The distribution of existing heat flow measurements is highly heterogeneous both on continents and on ocean basins [e.g., Pollack et al., 1993, Figure 1], so the recent upward correction is not surprising.

[14] Estimating the composition of the depleted mantle suffers from more compounded difficulties. At shallow depths, the depleted mantle usually differentiates into oceanic crust and depleted mantle lithosphere by melting, so the composition of the (unmelted) depleted mantle must be indirectly estimated based on these chemically differentiated phases. The oceanic crust itself is layered because of chemical fractionation, and only the very top portion (i.e., MORB) is usually accessible for sampling. So the estimate based on MORB has to take into account the effects of chemical fractionation as well as mantle melting. This indirectness may be compensated, however, by the fact that the oceanic crust is the product of mantle melting over the depths of several tens of kilometers. MORB samples thus have a potential to effectively probe the shallow upper mantle, with its small-scale heterogeneities automatically homogenized by melting. Perhaps a more serious issue is the considerable regional variations observed in the trace element chemistry of MORB, which tends to show the so-called “enriched” signature (higher concentrations of more incompatible elements) near hot spots such as Iceland [e.g., Schilling et al., 1983]. Samples from these anomalous regions have traditionally been disregarded when discussing the depleted mantle, because hot spots have often been considered to originate in a deep geochemical reservoir that is different from the depleted mantle. One may also be tempted to apply this type of screening to limit data variability and to estimate an average with small uncertainty. The origin of hot spot magmatism is, however, still debated (8), and this screening may not be warranted. The concentration of heat-producing elements in the depleted mantle has been estimated to be 8 ppb U, 16 ppb Th, and 100 ppm K by Jochum et al. [1983] and 4.7 ± 1.4 ppb U, 13.7 ± 4.1 ppb Th, and 60 ± 17 ppm K by Salters and Stracke [2004], but these estimates are biased to the depleted end-member of MORB samples because of this screening. An alternative approach by Workman and Hart [2005] based on the trace element systematics of abyssal peridotites (exhumed pieces of DML) may be able to avoid this type of bias, but it suffers from different kinds of indirectness; their approach depends on the BSE composition as well as assumed isotopic evolution due to continental growth. The 30% error assigned by Salters and Stracke [2004] translates to the heat generation of 2.9–5.5 TW, but it would be corrected upward if data screening is loosened [Langmuir et al., 2005]. At the same time, expanded data would certainly exhibit much more variability. A careful uncertainty analysis has not been published, but it would be vital to derive realistic bounds on heat production. An estimate without uncertainty is of little value when scrutinizing global mass balance.

[15] Finally, the composition of the bulk silicate Earth represents what we expect for the bulk silicate Earth. Compared to our estimates for the continental crust and the depleted mantle, therefore, it is most theoretical, depending critically on the validity of theoretical assumptions employed. A number of geochemists have estimated the BSE composition with different sets of assumptions (see a review by Lyubetskaya and Korenaga [2007a]), and the so-called pyrolite approach [e.g., McDonough and Sun, 1995] appears to be most robust as it requires only the following two assumptions: (1) the BSE should be found somewhere along the compositional trend displayed by mantle rocks, and (2) the ratios of refractory lithophile elements (RLE) in the BSE are chondritic. If either of these assumptions is invalid, it is impossible to estimate the BSE composition (not only by the pyrolite approach but also by other existing methods), so its uncertainty would be essentially unbounded. The first assumption implicitly requires whole mantle convection. Our mantle samples such as mantle xenoliths and massif peridotites are all from the upper mantle, so theoretically the composition of the lower mantle could be different from that of the upper mantle, and if so, the RLE ratios could not be imposed on the composition trend of the upper mantle samples because these cosmochemical constraints are valid only when considering BSE as a whole. For the estimated BSE composition to be justified, therefore, global mass balance must be self-consistent. This issue is explained in more detail in 4. The second assumption is reasonable because refractory elements have high condensation temperatures, so they are unlikely to have been fractionated from each other during planetary formation processes, and lithophile elements do not enter the core. The RLE ratios are the most stable ratios among various classes of chondrites, and, indeed, the isotope data of terrestrial samples require that at least Sm/Nd, Lu/Hf, and Th/U (these are all RLEs) should be within a few percent of the chondritic value [e.g., Lyubetskaya and Korenaga, 2007a]. With these assumptions, we can quantify the uncertainty of the BSE composition by taking into account the uncertainty of observed compositional trend as well as that of chondritic RLE ratios, and Lyubetskaya and Korenaga [2007a] obtained ∼20% error for the BSE abundance of refractory lithophile elements. This uncertainty represents the influence of scatter in the upper mantle compositional trend under the imposed chondritic constraints. Lyubetskaya and Korenaga also found that the previous estimate of 20 ppb U and 80 ppb Th should be corrected slightly downward, and the new estimate with one standard deviation is 17 ± 3 ppb U and 63 ± 11 ppb Th. As K is not a RLE, it must be derived with additional constraints such as bulk Earth K/U and K/La, which further amplify uncertainty. The revised estimate is 190 ± 40 ppm K. The likely range of BSE heat production is then 13–19 TW.

[16] Note that these uncertainties in heat production do not correlate to each other because the compositions of these silicate reservoirs are estimated independently of each other. Using one standard deviation therefore, the range of combined CC and DM heat production is 7.9–15.5 TW, which overlaps with the range of BSE heat production (note also that the upper bound of 15.5 TW is likely to increase if the DM composition is estimated with more even data sampling). The missing heat source paradox is thus unproven and probably does not exist; it could well be an artifact caused by neglecting unavoidable uncertainties in geochemical inference. Given the uncertainties, it is still possible to postulate a hidden reservoir (e.g., by assuming 5 TW for CC, 3 TW for DM, and 19 TW for BSE), but it is no longer required.

2.1.3. Importance of Self-Consistent Mass Balance

[17] The self-consistent mass balance is vital to justify the estimated BSE composition because, as noted in 3, whole mantle convection has to be assumed when deriving the composition. In other words, if mass balance is not self-consistent, the BSE composition is undefined. Interestingly, however, the notion of layered-mantle convection has coexisted with the BSE composition estimated by assuming whole mantle convection. A common layered-mantle model advocated by geochemists is illustrated in Figure 3b; the upper layer is the depleted mantle as observed through mid-ocean ridge magmatism, while the lower layer is identified as the primitive mantle [e.g., Jacobsen and Wasserburg, 1979; DePaolo, 1980; Allègre, 2002]. Note that the PM composition is identical to the BSE composition, both in major and trace elements, and is very similar to the DM composition in terms of major elements. The volume of the upper layer is determined so that the average composition of these three reservoirs (CC, DM, and PM) is equal to the BSE composition. This is, of course, one possibility, but the average composition can be (vastly) different from the estimated BSE composition both in major and trace elements (Figure 3c). In fact, if a layered-mantle model has to be true, the major element composition of the lower layer may be required to be different from the PM composition. In the classical layered-mantle model, the internal boundary coincides with the 660-km discontinuity, and it was speculated that the endothermic phase change expected at this discontinuity may be able to sustain the layering [e.g., Christensen and Yuen, 1984]. Later studies on mantle convection, however, suggested that it was difficult to maintain the layered state with the phase change alone [e.g., Tackley et al., 1993; Zhong and Gurnis, 1994], and this idea became obsolete particularly after the internal boundary was pushed down the 660-km discontinuity. To achieve a layered mantle in a physically plausible manner, therefore, a different major element composition (e.g., enriched in iron) seems necessary to make the lower layer intrinsically denser than the upper layer [e.g., Kellogg et al., 1999; Davaille, 1999].

[18] A failure to close global mass balance therefore presents a serious challenge to our understanding of the bulk Earth chemistry. It does not only indicate the presence of a hidden reservoir, but it also makes the BSE composition (and thus heat production) unconstrained (Figure 3c). The BSE composition is meaningful only with whole mantle convection (Figure 3a), and it cannot be used to estimate the composition of a hidden reservoir.

2.1.4. Missing Argon Paradox and More

[19] As discussed in 3, the missing heat source paradox does not seem to be as compelling as commonly described. Given the uncertainty involved, however, the argument could still be used to support layered-mantle convection if any of the other mass balance arguments is more robust, and it would be difficult for geochemistry to constrain the amount of internal heat production in the bulk of the mantle. Other arguments involve the mass balance of noble gases such as argon and xenon [e.g., Allègre et al., 1996; Porcelli and Turekian, 2003] and lithophile elements [e.g., Helffrich and Wood, 2001; Albarède and van der Hilst, 2002], but as recently reviewed by Lyubetskaya and Korenaga [2007b], on the basis of the revised BSE composition, these mass balances do not contradict whole mantle convection if uncertainty is properly taken into account.

2.1.5. Isotopic Mass Balance and Two Versions of Chondrite Assumptions

[20] The final kinds of mass balance to be discussed are the ones that involve isotopic ratios such as ¹⁴³Nd/¹⁴⁴Nd, ¹⁷⁶Hf/¹⁷⁷Hf, and, more recently, ¹⁴²Nd/¹⁴⁴Nd. Though similar to other mass balances discussed so far, isotopic mass balance arguments involve a much more strict definition for the chondrite assumption, warranting a separate treatment.

[21] The ¹⁴³Nd/¹⁴⁴Nd ratio, for example, reflects the long-term evolution of ¹⁴⁷Sm/¹⁴⁴Nd (or Sm/Nd) as ¹⁴³Nd is a decay product of ¹⁴⁷Sm. Chondritic meteorites exhibit a tight distribution of ¹⁴⁷Sm/¹⁴⁴Nd (0.196 ± 0.002), which corresponds to ¹⁴³Nd/¹⁴⁴Nd (0.51263 ± 0.00005) [Patchett et al., 2004, Figure 1] (Note that the ∼1% deviation in ¹⁴⁷Sm/¹⁴⁴Nd reduces to only ∼0.01% deviation in ¹⁴³Nd/¹⁴⁴Nd because of the long half-life involved.) The ¹⁴³Nd/¹⁴⁴Nd ratio of terrestrial samples is usually compared to this well-defined (present day) chondritic value or its value in the past, and the difference is expressed in the ε¹⁴³Nd notation defined as

where t is the age of a sample and CHUR refers to a chondritic uniform reservoir. A common argument is that the BSE ε¹⁴³Nd would be ∼5 if CC contains 20 ppm Nd with ε¹⁴³Nd of −16 and DM contains 0.7 ppm Nd with ε¹⁴³Nd of 9 [e.g., Lassiter, 2004]. A major source of uncertainty is the depleted mantle, which stores more than 80% of total Nd in this mass balance. Though both ε¹⁴³Nd and Nd content of MORB vary from sample to sample by ∼30%, they tend to be negatively correlated (i.e., a more enriched sample has a lower ε¹⁴³Nd) [e.g., Hofmann, 2003, Figure 14], so it is difficult to lower the BSE ε¹⁴³Nd below ∼3 even if we take into account all of the likely uncertainties (though a more complete statistical study of MORB chemistry may indicate otherwise in the future). Because chondritic samples show ε¹⁴³Nd of ±1 as mentioned above, this relatively high BSE value seems to require a hidden reservoir with low Sm/Nd (i.e., enriched in incompatible elements because Nd is more incompatible than Sm).

[22] Note that the chondrite assumption for BSE is more strict here than in the case of estimating the BSE composition. Requiring the BSE ε¹⁴³Nd to be within ±1 is equivalent to assuming that the ratio Sm/Nd is within 1% deviation from the chondrite average, but the Sm-Nd isotopic evolution recorded by terrestrial samples does not place such a tight constraint [Lyubetskaya and Korenaga, 2007a, 2007b]. Up to a few percent deviation is acceptable because the earliest record for the depleted mantle at ∼4 Ga ago shows ε¹⁴³Nd of as high as ∼4 [e.g., Bennett, 2003]. Thus, requiring the RLE ratios such as Sm/Nd to be within a few percent deviation from the chondritic value is justified by existing data, whereas imposing the ±1% deviation is not. It is important to distinguish between these weak and strong versions of the chondrite assumptions. When the strong version is not satisfied by a certain mass balance, it is often called “nonchondritic,” and this might imply that refusing the presence of a hidden reservoir violates cosmochemical considerations. We do not understand well, however, how efficiently elements and isotopes were homogenized in the presolar nebula and how much they were later modified by planetary accretion processes [e.g., Halliday and Porcelli, 2001].

[23] As already mentioned, even ¹⁴⁷Sm/¹⁴⁴Nd in chondrites has as much as ±1% variation (and ¹⁷⁶Lu/¹⁷⁷Hf has more than ±5% variation [Patchett et al., 2004]). Furthermore, chondrites are believed to originate mostly from the asteroid belt [e.g., Scott and Krot, 2003], so even though the number of chondrite samples is large, we are essentially looking at the highly restricted portion of the solar system; the mass of the asteroid belt is less than 0.1% of that of Earth. The stable statistics of chondrite compositions may be just due to more or less repetitive sampling. Earth is the best sampled planet, and yet we cannot measure its bulk isotope characteristics. We have only a few dozen meteorites from Mars, and no sample for Venus and Mercury, so their bulk isotopic compositions are even more uncertain than Earth's. Note that the spectroscopic measurement of the photospheric composition of the Sun may have ∼25% uncertainty regarding Sm/Nd [Grevesse and Sauval, 1998] (assuming that reported errors are not correlated). A remarkable correlation between the photospheric composition and the chondrite composition exists for a number of elements, but this correlation in the logarithmic space does not impose a tight constraint on the ratios of trace elements. Some authors are well aware of this possibility of limited sampling by chondrites [e.g., Wasson and Kallemeyn, 1988; Drake and Righter, 2002], but the strong version of the chondrite assumption has long been a part of central dogma in geochemistry.

[24] Recently, Boyet and Carlson [2005] proposed a new kind of argument for a hidden reservoir based on the now extinct isotope ¹⁴⁶Sm, which decayed into ¹⁴²Nd. So far, measured terrestrial samples show systematically higher ¹⁴²Nd/¹⁴⁴Nd than chondrites by ∼0.2 ± 0.14 ε unit. To have the BSE ε¹⁴²Nd be (strictly) chondritic, therefore, there must be an unobserved reservoir, whose ε¹⁴²Nd is more negative than chondritic. Also, because of the short half-life of ¹⁴⁶Sm, this reservoir has to have formed within the first 100 Ma of Earth's history. This presumed reservoir does not seem to be sampled even by hot spot magmatism, so if its volume is substantial, it must somehow be able to avoid entrainment by the convecting mantle. Alternatively, the reservoir could be trivially small or even nonexistent if the BSE ε¹⁴²Nd (or Sm/Nd) does not have to exactly equal the chondritic average. With igneous differentiation, Sm/Nd in the outer portions of planetesimal bodies would be low, and impact erosion could preferentially remove major portions of the outer silicate portions of accreting planetesimals, resulting in a slightly elevated Sm/Nd ratio in the growing Earth [Halliday, 2003]. As Lyubetskaya and Korenaga [2007b] argued, both ε¹⁴³Nd and ε¹⁴²Nd arguments rely heavily on the strong version of the chondrite assumption. Andreasen and Sharma [2006] showed that ordinary chondrites have systematically higher ε¹⁴²Nd than carbonaceous chondrites, and the bulk Earth may be located somewhere along this trend. To summarize, given the resolution of available observational constraints on the BSE composition, existing mass balance arguments, whether elemental or isotopic, do not necessarily contradict the notion of whole mantle convection.

2.2.1. Origins of Hot Spots

[27] Large-scale mantle circulation associated with plate tectonics does not predict the existence of hot spots, so some kind of small-scale anomalies are necessary to induce mantle melting. The plume hypothesis provides thermal anomalies in the form of convective instabilities rising from the bottom boundary layer. Though such instabilities are likely to exist in the cooling Earth, whether plumes can explain all hot spots is another issue. At least for some hot spots, in fact, the chemistry of hot spot magma is clearly inconsistent with the melting of hot mantle, suggesting that chemically anomalous mantle may play an important role [e.g., Korenaga and Kelemen, 2000]. In recent years, debates over the existence of mantle plumes have been particularly intense [e.g., Anderson, 2000; DePaolo and Manga, 2003; Foulger and Natland, 2003; Sleep, 2003; Campbell, 2005; Foulger et al., 2005b; McNutt, 2006] but sometimes appear to be largely philosophical and semantic. What is frustrating is that observations tend to be ambiguous enough to allow conflicting interpretations to coexist and that our understanding of mantle dynamics is too incomplete in a number of aspects to allow robust predictions. It is perhaps true that the notion of mantle plumes has been overused in the past, but it would also be extreme to reject it entirely. Some hot spots such as Iceland and Galapagos may require some special mechanism other than a mantle plume [e.g., Korenaga, 2004; Sallares et al., 2005], but others such as Hawaii may simply be caused by plumes. It might be aesthetically pleasing to have just one theory to explain all hot spots, but existing observations indicate that the reality is unfortunately rather complex [e.g., Courtillot et al., 2003].

[28] The Iceland hot spot, for example, is a classical example of a ridge-centered hot spot, but it does not fit well into the standard theory of a mantle plume. A number of seismological studies show that Iceland is underlain by a low-velocity anomaly at least in the shallow upper mantle [e.g., Wolfe et al., 1997; Foulger et al., 2001], but global tomography indicates that the anomaly does not extend into the lower mantle [e.g., Ritsema et al., 1999; Montelli et al., 2004]. (Note that some studies have presented a tomographic image for a deep Icelandic plume [Bijwaard and Spakman, 1999; Zhao, 2001], but it appears to be an artifact of the use of an oversaturated color scale [Foulger and Pearson, 2001].) Moreover, the geochemistry of Icelandic lavas suggests that the excess magmatism may be caused (at least partially) by an unusually fertile source mantle because of recycled oceanic crust [Korenaga and Kelemen, 2000; Foulger et al., 2005a; Sobolev et al., 2007]. Major element heterogeneity of this sort is not unexpected in the convecting mantle. At mid-ocean ridges, plate tectonics continuously differentiates the depleted mantle into the oceanic crust and the depleted mantle lithosphere, and these two phases are remixed sometime after subduction. However, subducted oceanic crust could segregate from the mantle lithosphere and advect independently [e.g., Karato, 1997; Korenaga, 2004; Lee and Chen, 2007], so premature remixing may lead to regional enrichment or depletion of the crustal component. The compositional variability observed along the mid-ocean ridge system (including anomalous regions such as Iceland) may simply reflect how remixing actually takes place in Earth's mantle.

[29] Note that fertile mantle is usually denser than normal, and it is expected to sink unless the chemical buoyancy is compensated by the thermal buoyancy [e.g., O'Hara, 1975; Cordery et al., 1997]. Even without thermal anomalies, however, dense fertile mantle may be brought up to the surface if, for example, it is located near divergence plate boundaries [e.g., Korenaga, 2005b]. Obviously, such a mechanism does not apply to midplate hot spots such as Hawaii. Though other alternatives have been proposed for midplate hot spots [e.g., Smith, 2003], it should be kept in mind that plume-like upwelling is hard to avoid in the cooling Earth, so at least some hot spots are expected to originate in deep mantle plumes. As briefly mentioned in 1 and discussed further in 12, the heat production of the convecting mantle is far from being sufficient to explain the observed heat flux, so Earth's mantle must have been cooling, which is also consistent with petrological data [e.g., Abbott et al., 1994; Grove and Parman, 2004]. The core must then be cooling, and so we have finite core heat flux. This means that the mantle cannot be purely internally heated, and it should also be heated from below. Core cooling is also fundamental to maintain the geodynamo. The hot thermal boundary layer at the core-mantle boundary is therefore a corollary of the present-day heat budget, unless all of core heat flux is consumed to just warm up cold slabs arriving at the core-mantle boundary. The convective instability of this boundary layer is only one step further in this line of logic. With temperature-dependent viscosity, the instability should evolve into a large plume head followed by a narrow conduit [e.g., Whitehead and Luther, 1975; Griffiths and Campbell, 1990].

[30] The form of convective instability, however, depends strongly on the rheology of the lowermost mantle and may be different from what the standard plume theory assumes [Korenaga, 2005a; Čadek and Fleitout, 2006]. Whereas the rheology of the upper mantle is now understood reasonably well [e.g., Karato and Wu, 1993; Hirth and Kohlstedt, 2003; Korenaga and Karato, 2008], that of the lower mantle is still largely unknown because of experimental difficulties. Only a broad agreement has been achieved for spatially averaged viscosity on the basis of geodynamical estimates [e.g., King, 1995; Forte and Mitrovica, 2001], and any other details such as temperature and pressure dependency are yet to be determined. The dynamics of the core-mantle boundary region may be further complicated by the likely presence of chemical heterogeneities [e.g., Farnetani, 1997; Montague and Kellogg, 2000; Jellinek and Manga, 2002], the postperovskite phase transition [e.g., Murakami et al., 2004; Nakagawa and Tackley, 2004], and drastic changes in thermal properties [e.g., Badro et al., 2004; Lin et al., 2005; Matyska and Yuen, 2005], so even if we limit ourselves to the dynamics of thermal plumes to a static background, our understanding is still immature. Note that this directly translates to the uncertainty of core-mantle interaction, i.e., how efficiently the mantle could cool the core (19).

2.2.2. Mantle Viscosity, Mantle Mixing, and Sampling by Magmatism

[31] Numerical studies on mantle mixing indicate that it is difficult to preserve a large-scale geochemical reservoir over Earth's history in the convecting mantle by anything other than density stratification, such as phase transitions, higher viscosity in the deep mantle, and/or depth-dependent thermal properties [e.g., van Keken and Ballentine, 1998, 1999; Hunt and Kellogg, 2001; Naliboff and Kellogg, 2007]. Though the presence of multiple geochemical reservoirs with different isotopic characteristics is a well-established observational fact, most notably supported by the studies of noble gas isotopes [e.g., Kurz et al., 1982; M. Honda et al., 1993; Tolstikhin and Hofmann, 2005], this observation by itself does not necessarily demand a large-scale reservoir, as noted in 8, and the reevaluation of global mass balance arguments casts a doubt on the existence of a massive reservoir to begin with (2). The difficulty in preserving a large-scale reservoir thus does not present any challenge to reconciling mantle dynamics and geochemistry.

[32] In contrast, small-scale reservoirs, such as distributed blobs in the so-called “marble cake mantle” or “mantle plum pudding” [e.g., Allègre and Turcotte, 1986; Phipps Morgan and Morgan, 1999], are expected to be better preserved [e.g., Christensen and Hofmann, 1994; Davies, 2002; Xie and Tackley, 2004], and thus whole mantle convection with small-scale chemical heterogeneities (Figure 2b) may be a promising concept both geophysically and geochemically. The fate of small-scale heterogeneities in the convecting mantle and their relation to geochemical observations are, however, difficult issues to be fully explored in future. A few key issues are highlighted in the following.

[33] First of all, “small-scale” heterogeneities refer here to anything but a massive reservoir that occupies a significant portion of the mantle, so their spatial scales are only vaguely defined. Theoretical studies on mantle mixing suggest that convective mixing, through repeated stretching and folding, could easily homogenize heterogeneities down to submeter scales within a billion years or so [e.g., Olson et al., 1984; Hoffman and McKenzie, 1985; Ferrachat and Ricard, 1998], but note that such studies focus on the destruction of existing heterogeneities by isoviscous convection. Though isotopic heterogeneities are most easily detected by geochemistry and thus most frequently discussed, they are likely to be accompanied by variations in trace element and major element concentrations (e.g., subducted oceanic crust and depleted mantle lithosphere), which, in turn, should accompany viscosity variations because viscosity is a strong function of lithology as well as trace elements such as hydrogen [e.g., Karato, 1986; Hirth and Kohlstedt, 1996; Karato, 1997]. Convective mixing in the presence of such composition-dependent viscosity is expected to be inefficient [Manga, 1996]. Moreover, as noted in 8, chemical heterogeneities are continually generated by plate tectonics, which must adversely affect homogenization. Field observations indeed indicate that chemical heterogeneities exist at almost all plausible scales: from centimeters to more than a few hundreds kilometers [e.g., Zindler and Hart, 1986; Carlson, 1994]. A number of seismological studies also support the existence of chemical heterogeneities throughout the mantle [e.g., Hedlin et al., 1997; Helffrich and Wood, 2001]. These observations are not surprising because a chemically heterogeneous mantle is a rule, not an exception, in the presence of plate tectonics, but at the same time they demand that we tackle a highly complex problem. For example, we know that a wide range of spatial scales must be involved, but we do not know the relative significance of different scales [cf. Gurnis, 1986]. The nature of chemical heterogeneities could also vary from place to place [e.g., Hirschmann and Stolper, 1996]. The concept of chemically heterogeneous mantle thus necessarily involves a few poorly constrained parameters, so unfortunately it tends to be regarded as an ad hoc explanation or a “just so” story, especially if applied to regional problems.

[34] Another important issue regarding chemical heterogeneities is how they are sampled by geochemistry. Different lithologies are likely to have not only different viscosities but also different densities, and as mentioned in 8, more fertile mantle may be preferentially sinking and thus less likely to be sampled by melting, which usually requires upwelling. The upwelling of intrinsically dense heterogeneities requires special conditions such as thermal anomalies and plate tectonic stresses, and the interplay among these competing factors needs to be systematically explored. Also, how a composite mantle with different lithologies would melt is an interesting problem by itself, which is still poorly understood. Enriched crustal fragments embedded in a more depleted mantle matrix, for example, are expected to melt more easily than the matrix, and because of this differential melting behavior, a single composite mantle could potentially yield melts with different isotopic characteristics depending on whether its upwelling (and thus melting) is limited to the base of lithosphere (in case of midplate hot spots) or not (under mid-ocean ridges) [e.g., Sleep, 1984; Meibom and Anderson, 2003; Ito and Mahoney, 2005]. Geochemical evidence for such composite mantle has been slowly accumulated [e.g., Hauri, 1996; Takahashi et al., 1998; Korenaga and Kelemen, 2000; Sobolev et al., 2007], but inferring the composition of a source mantle and its potential temperature from lava samples is still a daunting task because predicting a lava composition from a given composite mantle is much more involved than the case for a homogeneous mantle. In a homogeneous mantle, for example, the melt phase can almost always migrate through the residual mantle because a porous flow network is expected to form at a very low porosity [Kohlstedt, 1992], but it is not clear whether the melting of a composite mantle behaves in a similar manner. Early melts in the enriched fragments may not escape into the surrounding matrix until the matrix itself starts to melt, and because the early melting extracts heat from the matrix [Sleep, 1984; Phipps Morgan, 2001], it could also delay the melting of the matrix. The consequence of polybaric melting such as chemical reaction between the residual mantle and the migrating melt phase is already a complex problem even for a homogeneous mantle [e.g., Kelemen et al., 1995; Spiegelman and Kelemen, 2003], and extending it to the case of a composite mantle may require more than a brute force approach given compositional diversity and various spatial scales involved.

[35] Finally, duration is an essential factor for mixing. Even in the case of simple isoviscous mixing, for example, homogenization would not be achieved if convection lasts for the period of only 100 Ma. Earth's mantle has been convecting probably for more than 4 Ga, but given the variable mantle properties, can it be considered sufficiently long? A good starting point is to estimate how quickly the mantle can be processed by mid-ocean ridge magmatism, which may be expressed as

where ρ_m is the density of shallow upper mantle (3300 kg m⁻³), R is the rate of plate construction at mid-ocean ridges, and D is the initial depth of mantle melting. With the current plate motion and mantle temperature, the processing rate is estimated to be ∼6.7 × 10¹⁴ kg a⁻¹ [Korenaga, 2006], so it would take ∼6 Ga to process the entire mantle. It is commonly assumed that mantle convection was more vigorous in the past, so this estimate based on the present-day processing rate may not be useful. In fact, some recent mixing calculations explicitly simulate more rapid convection in the past [e.g., Davies, 2002; Xie and Tackley, 2004]. Greater convective vigor in the past is, however, merely an assumption without geological evidence, and it is also known to cause a serious problem in the thermal history of Earth (12). A scaling law for plate tectonic convection suggests that mantle convection may have been less vigorous when the mantle was hotter in the past [Korenaga, 2003, 2006], and if so, the present-day processing rate could be a good approximation for the last 4 Ga or so (Figure 5f). This slow processing rate may be a good proxy for how quickly chemical heterogeneities are created and destroyed by terrestrial magmatism in general. In order to predict how efficiently convection could homogenize the mantle over the geological time, therefore, we need to have a good idea of how mantle convection has evolved over the last 4.6 Ga. Simulating realistic mantle mixing requires successfully simulating the long-term evolution of mantle convection at the same time, but doing so without violating observational constraints on the thermal budget is yet to be attempted.

[36] It should be clear now that we have just begun to understand how geochemical observations and mantle convection may be related. A deep, isolated reservoir has often been called for to preserve isotopic heterogeneities, but it is based mostly on the following two assumptions: (1) hot spot magmatism is caused by deep mantle plumes, and (2) mantle convection can quickly homogenize chemical heterogeneities. Both of them need to be carefully examined in the framework of plate tectonic convection. Chemically heterogeneous mantle, in terms of isotopes as well as trace and major elements, is an unavoidable consequence of mantle convection with plate tectonics, carrying a number of important implications for terrestrial magmatism and mantle dynamics. Whole mantle convection is conceptually a simple mechanism, but it can create complex phenomena, so its potential to reconcile geophysics and geochemistry should not be underestimated.

3.3.1. Variations on Layered-Mantle Models

[48] As mentioned in 4, the endothermic phase change expected at the 660-km discontinuity was once considered a promising mechanism to achieve the particular kind of layered mantle favored by geochemists (Figure 3b), and numerical studies on mantle convection suggested that layering would be possible if the Clapeyron slope is more negative than −6 MPa K⁻¹ [Christensen and Yuen, 1984, 1985]. Later, experimental studies indicated a less negative (but still significant) Clapeyron slope for this phase change (∼−2.8 MPa K⁻¹) [e.g., Ito and Takahashi, 1989], and this reinvigorated an interest in the role of the endothermic phase change in mantle dynamics [e.g., Machetel and Weber, 1991; Tackley et al., 1993; S. Honda et al., 1993; Solheim and Peltier, 1994; Yuen et al., 1994]. These studies pointed to the possibility of episodically layered convection, and in particular, Yuen et al. [1994] suggested that the degree of layering may increase with an increasing Rayleigh number, which implies that the mantle could have been layered when the vigor of convection was greater because of hotter and less viscous mantle. Motivated by this, Honda [1995] explored the implications this transition with parameterized convection models, but his conclusion is very similar to that of Christensen [1985a]; even with the transition, an unconventional heat flow scaling must be assumed to satisfy the Urey ratio constraint. Davies [1995] investigated the consequence of episodic layering on thermal history but with the present-day Urey ratio of >0.8. More recently Butler and Peltier [2002] developed a more elaborate parameterization for time-dependent layering, and they found that a gradual (not sudden) transition from layered to whole mantle convection could help to satisfy the Urey ratio constraint while keeping a conventional heat flow scaling.

[49] Recent experimental studies indicate, however, that the Clapeyron slope of the endothermic phase change is likely to be much less negative than previously thought (now estimated to be around −1.3 MPa K⁻¹) [Katsura et al., 2003; Fei et al., 2004]. Furthermore, the role of the endothermic phase change tends to be overemphasized in a majority of previous numerical studies by the use of unrealistic mantle viscosity (such as constant viscosity [e.g., Yuen et al., 1994]); convection is characterized by short wavelengths, lacking long-wavelength structures associated with plate tectonics. This is important because the degree of layering depends not only on the Clapeyron slope but also on the wavelength of convective flow [Zhong and Gurnis, 1994; Tackley, 1995]; large and stiff plates prefer whole mantle convection. On the basis of the chemical buoyancy of oceanic plates, Korenaga [2006] argued that plate dimensions in the past were unlikely to be very different from the present, implying a diminishing chance for time-dependent layering.

3.3.2. Variations on Whole Mantle Models

[50] As the continental crust today probably contains about 40% or more of total heat-producing elements in the bulk silicate Earth (3), its evolution in the past must have affected the thermal history of Earth to some degree. Also, continents are commonly thought to “insulate” Earth's mantle [e.g., Anderson, 1982; Gurnis, 1988], so they may reduce the degree of secular cooling. Spohn and Breuer [1993] investigated these two effects by extending the parameterized convection model of Stevenson et al. [1983] and suggested that the influence of continental growth would be minor, reducing the predicted value of the present-day Urey ratio by only ∼0.1. Various growth scenarios simulated by Spohn and Breuer [1993] are, however, all similar to the so-called “instantaneous” growth model of Armstrong [1981], in which the present-day continental volume is achieved in the early Archean and is then maintained by a perfect balance between the creation of new crust and the recycling of old crust. Grigné and Labrosse [2001] revisited this issue with a more flexible parametrization of continental growth, but the other aspects of their modeling are rather tightly fixed. In particular, their modeling strategy was built solely on the heat flow scaling law developed for an isoviscous fluid [Sotin and Labrosse, 1999], and other possibilities were not explored. The limited parameter search impairs the robustness of their conclusions about a likely scenario for continental growth. They also modeled continents as perfect insulators, but whether continents actually insulate the mantle in terms of global heat loss has recently been called into question [Lenardic et al., 2005]. Given the uncertainty involved, thermal history calculations are unlikely to place a tight constraint on continental growth [e.g., Lenardic, 2006]. Nonetheless, incorporating a plausible growth history is important to improve the accuracy of a predicted thermal history because it definitely affects the heat production of the convecting mantle. The history of continental growth and recycling is still highly controversial among geochemists [e.g., Collerson and Kamber, 1999; Campbell, 2003; Harrison et al., 2005; Hawkesworth and Kemp, 2006], but if a consensus is ever achieved, thermal evolution modeling may become useful for inferring Archean and Hadean dynamics [e.g., Korenaga, 2006].

[51] The strongly temperature-dependent nature of the heat flow scaling law (equation (9) with m ∼ 10) originates in high activation energies measured for the deformation mechanisms of olivine, which is stable only in the upper mantle (i.e., above the 410-km discontinuity). Whereas the average viscosity of deeper mantle under present-day conditions can be inferred from geophysical observations such as the geoid and postglacial rebound [e.g., Hager et al., 1985; Mitrovica and Forte, 1997], its temperature dependency must be estimated through experimental studies, which are still difficult to conduct at lower mantle conditions [e.g., Karato et al., 1995]. Also, we can only speculate at this point how different deformation mechanisms, such as diffusion and dislocation creep, might compete in the lower mantle [Karato, 1998], and van den Berg and Yuen [1998] suggested that, under certain conditions, heat flow scaling for the lower mantle may resemble that of Christensen [1985a]. Alternatively, Solomatov [1996] suggested that if the lower mantle is mainly deformed by diffusion creep, which is known to be very sensitive to grain size, the kinetics of grain growth may reverse the temperature dependency of lower mantle viscosity; that is, hotter mantle may become stiffer because grains would grow faster at higher temperatures. Solomatov [2001] pursued this idea with parameterized convection models, but his conclusions were unavoidably equivocal because the kinetics of grain growth itself is still under debate [e.g., Yamazaki et al., 1996; Solomatov et al., 2002] and how grain growth should interact with convective deformation is yet to be understood [e.g., Hall and Parmentier, 2003; Bercovici and Ricard, 2005]. Recently, however, Korenaga [2005a] pointed out that the present-day plume dynamics may provide an observational support for this grain-size-dependent rheology.

[52] Yet another mechanism to delay secular cooling was proposed by van den Berg and Yuen [2002] based on variable thermal conductivity. As temperature increases, the thermal conductivity of mantle minerals first decreases rapidly and then might increase gradually owing to the competing effect of conduction and radiation processes [Hofmeister, 1999], and this temperature dependency can result in a low-conductivity zone in the shallow upper mantle (though this dependency is controversial [cf. Shankland et al., 1979]). van den Berg et al. [2004] extended their earlier study with constant viscosity to temperature- and pressure-dependent viscosity and simulated the thermal evolution of Earth. In their studies, the case of variable conductivity is compared with that of constant conductivity to demonstrate less efficient cooling for the former, but this is simply because the constant conductivity they used corresponds to the surface value of variable conductivity (i.e., the highest conductivity). Note that the heat flow scaling of the kind of equation (9) does not assume a particular value of thermal conductivity because the absolute value of its prediction is normalized with the present-day heat flux. The relative variation in predicted heat flow (as a function of internal temperature) is what matters most and is governed largely by the temperature dependence of mantle viscosity, compared to which the effect of variable conductivity is negligible. A moderate degree of secular cooling observed in the simulation of van den Berg et al. [2004] is because mantle viscosity is assumed to be weakly temperature-dependent in their calculation.

[53] Recently, by combining several geological and geochemical arguments, Silver and Behn [2008] proposed that the efficiency of plate tectonics may have experienced large fluctuations (as much as an order of magnitude) on a billion-year timescale and suggested that such intermittent plate tectonics may help to avoid thermal catastrophe even with the conventional heat flow scaling. Their heat flow parameterization, however, fails to account for the significance of slow thermal diffusion and the role of convective instability in the limit of stagnant lid convection, and their solution for thermal evolution can be shown to be physically unrealistic. They suggested an order-of-magnitude reduction in heat flux at ∼1 Ga ago, which should be accompanied by an order-of-magnitude increase in the average thickness of the top thermal boundary layer (or oceanic lithosphere); lithosphere would then occupy the entire upper mantle. Given our understanding of upper mantle rheology, convective instability would certainly prevent such thick lithosphere, and even if the mantle is rigid and allows the indefinite growth of the boundary layer, it would take as long as about 2 Ga and not be instantaneous as assumed by Silver and Behn [2008]. Intermittent plate tectonics by itself, even if true, does not help to solve the long-standing paradox in Earth's thermal evolution, and available geological records for the Phanerozoic eustasy [Algeo and Seslavinsky, 1995; Miller et al., 2005] indicate that its influence on surface heat flux should be small [cf. Korenaga, 2007].

3.3.3. Plate Tectonics and Whole Mantle Convection

[54] Efforts to reconstruct thermal evolution by starting from establishing heat flow scaling appropriate for Earth's mantle may be called a physical approach; its underlying principle is that a better understanding of the physics of mantle convection should also help to build a reasonable thermal history. In contrast, one could also infer an acceptable form of heat flow scaling simply by imposing observational constraints on thermal evolution. For example, because thermal catastrophe results from the strongly temperature-dependent heat flow scaling, it is obvious that such an unrealistic history can be avoided if convective heat flow does not depend much on internal temperature. A similar conclusion can also be obtained using the argon budget [Sleep, 1979; Tajika and Matui, 1993]. Conclusions drawn from this empirical approach are essentially the same as those of Christensen [1985a], and the only difference between these two approaches is whether one strives to understand the physics behind heat flow scaling. Labrosse and Jaupart [2007] have recently revived the empirical approach by noting that there is no theoretical framework currently available to predict an area-age relation for seafloor, which is essential to obtain an accurate estimate of convective heat flux. The empirical approach, however, has no predicting power because a physical mechanism is left unresolved. Convection in Earth's mantle is indeed a complex phenomenon, and it may be tempting to get discouraged, but it would be unwise to settle for an ad hoc solution by bypassing the important physics.

[55] The operation of plate tectonics is one of the unique features of Earth, but it is not fully understood why this special kind of mantle convection can take place to begin with [e.g., Tackley, 2000a; Schubert et al., 2001; Bercovici, 2003]. With the strongly temperature-dependent viscosity of Earth's mantle, the most likely style of thermal convection should be stagnant lid convection [Solomatov, 1995], in which convection takes place only beneath a rigid spherical shell. On Earth, however, the surface is broken into a number of rigid plates, which can recycle back to deep mantle by subduction. Besides temperature-dependent viscosity, therefore, there must be some additional mechanism to weaken the otherwise stiff surface, and how to simulate plate tectonic convection by introducing some kind of self-sustaining weakening has been intensively studied particularly since the mid-1990s [e.g., Bercovici, 1996; Moresi and Solomatov, 1998; Trompert and Hansen, 1998; Tackley, 2000b; Richards et al., 2001; Auth et al., 2003; Gurnis et al., 2004; Stein et al., 2004; Bercovici and Ricard, 2005].

[56] Along with these studies, our understanding of the likely energy balance in plate tectonic convection has been gradually improved. Solomatov [1995] first noted that the deformation of a cold and strong boundary layer (i.e., the bending of a subducting plate) must accompany a significant fraction (∼50%) of total energy dissipation, and this was investigated in further detail by Conrad and Hager [1999a, 2001]. On the basis of the notion that plate bending may be the bottleneck of plate tectonics, Conrad and Hager [1999b] explored its implications for heat flow scaling and suggested that a temperature-independent scaling may be possible. Their scaling, however, depends on assumed plate thickness, to which bending dissipation is very sensitive. The maximum thickness of oceanic plates is usually considered to be ∼100 km [e.g., Parsons and Sclater, 1977; Stein and Stein, 1992; McKenzie et al., 2005], but this is merely the present-day value. To what thickness oceanic plates can grow under different conditions was not understood well until recently [e.g., Korenaga and Jordan, 2002, 2003; Huang et al., 2003; Zaranek and Parmentier, 2004]. On the basis of this better understanding of plate thickness variation, Korenaga [2003] revisited the scaling of Conrad and Hager [1999b] and found that bending dissipation alone would not reduce the temperature dependence because plates should get thinner for a hotter mantle. Plate tectonics is, however, accompanied by mantle melting beneath divergent boundaries, and this melting strengthens the depleted mantle lithosphere because of dehydration stiffening [e.g., Hirth and Kohlstedt, 1996]. As a hotter mantle starts to melt at a greater depth, this results in a thicker depleted lithosphere [e.g., McKenzie and Bickle, 1988]; therefore, it is possible to expect a thicker plate and thus lower heat flux for a hotter mantle, reversing the sense of conventional temperature dependency [Korenaga, 2003]. The implications of this scaling for the thermal evolution of Earth have been explored extensively, and it is shown to be compatible with a range of existing geological, geophysical, and geochemical observations [Korenaga, 2006; Lyubetskaya and Korenaga, 2007b].

[57] Reconstructed thermal histories according to the formulation of Korenaga [2006] are shown in Figure 5 for the convective Urey ratio of 0.23 ± 0.15, with and without continental growth. To be consistent with the uncertainty analysis in 12, the convective heat flux and the continental heat production at present are set to 36 and 10 TW, respectively, for the Urey ratio of 0.08 and to 37 and 5 TW for the Urey ratio of 0.38. Linear continental growth starting at 4 Ga ago and continuing to the present is assumed. There is no reason to prefer this particular growth history over other possibilities; it is used here for illustrative purposes only. As in Figure 4, the thermal histories are integrated backward in time. The uncertainty of the Urey ratio (0.08–0.38) results in an ∼200 K difference at 3 Ga ago (Figure 5a). Though this temperature ambiguity is certainly nontrivial, one might also see it as surprisingly small given the large uncertainty of the Urey ratio (see/compare Figure 4a). This is because of a negative feedback originating in the heat flow scaling of Korenaga [2003], which also enables stable predictions for heat flow, plate velocity, and mantle processing rate in the past (Figures 5b, 5e, and 5f). Note that the hotter and stiffer lower mantle proposed by Solomatov [1996, 2001] would reinforce this negative feedback. The effects of continental growth can be seen as minor, as already indicated by previous studies [e.g., Spohn and Breuer, 1993].

[58] Recently, Labrosse and Jaupart [2007] suggested that almost all of previous heat flow scaling laws obtained by the physical approach do not apply to Earth's mantle because they are based on a single convection cell model and thus fail to capture the observed area-age distribution of seafloor. This argument is, however, not consistent with the nature of heat flow “scaling” in parameterized convection. Heat flow scaling is usually normalized such that it correctly yields the present-day heat flow with the present-day internal temperature, and what the scaling tries to predict is how this heat flux would change when the internal temperature is varied. In other words, the present-day area-age distribution is implicitly adopted as a reference convective planform, and its sensitivity to temperature perturbations is assessed by a simplified convection model. Whether this approach is too simplified to be useful is debatable. The conventional scaling of Schubert et al. [1980] and Davies [1980], for example, usually predicts very high heat flux and rapid plate tectonics in the past (Figures 4b and 4e), implying corresponding flow geometry that may be vastly different from today. In this case, deviations from the reference may be too large for a simple scaling to be useful. Conversely, the scaling of Korenaga [2006] indicates that the style of plate tectonics may not have changed much in the past (Figures 4b and 4e), which provides a posteriori justification for this scaling.

[59] The normalization with the present-day heat flux implies that the derived scaling is useful only for plate tectonic convection, so the operation of plate tectonics has to be assumed when reconstructing a thermal history. Plate tectonics was probably present at least back to the late Archean (∼3 Ga ago) [e.g., Hoffman, 1997], but the emergence of plate tectonics in Earth's history is not understood well. The very early Earth is likely to have been governed by a different mode of mantle convection [e.g., Sleep, 2000; Stevenson, 2003]. It is thus preferred to integrate a thermal history backward from the present. As emphasized by Korenaga [2006], a common forward integration is not recommended because it has to assume the same physics of convection throughout Earth's history; a more recent part of the history is unavoidably affected by this strong assumption of Hadean and Archean dynamics.

[60] The backward reconstruction of a thermal history, as shown in Figure 5, is in a sense promising because it may provide important boundary conditions on mantle dynamics in the early Archean and the Hadean. The negative feedback suggested by the recent heat flow scaling already makes robust predictions for heat flow and plate velocity in the past. Through its unifying framework, thermal evolution modeling assimilates a variety of geophysical and geochemical observations and converts them to quantitative constraints on the distant past. A better understanding of the physics of plate tectonics, a more accurate estimate on the heat production of the convecting mantle, and a less ambiguous history of continental growth will all help to sharpen this theoretical inference.