Lithologic Controls on Silicate Weathering Regimes of Temperate Planets

The atmospheres of temperate planets may be regulated by geochemical cycles. Silicate weathering provides essential negative feedback to the carbonate-silicate cycle (carbon cycle) to maintain temperate climates on Earth and possibly on Earth-sized temperate exoplanets. The intensity of weathering is normally attributed to the kinetics of weathering reactions of individual minerals. The implementation of a transport-controlled weathering model shows that when the CO$_2$ volume mixing ratio decreases or surface temperature increases, equilibrium thermodynamics rather than kinetics exerts a strong control on weathering. Modeling the weathering of all minerals in a given rock instead of individual minerals is crucial. The transition between kinetically- and thermodynamically-limited regimes of weathering is strongly sensitive to rock lithology. Application of this model to Earth suggests that global mean continental granite and seafloor basalt weathering rates are likely limited by the supply of fresh rocks, yet regional weathering rates can be influenced by both kinetics and thermodynamics. Consideration of total dissolved inorganic carbon as a proxy for weathering results in another CO$_2$ drawdown regime: CO$_2$ dissolution, where aqueous bicarbonate and carbonate ions produced by rock weathering are lower in concentration than aqueous CO$_2$. Upper limits to weathering as a function of lithology are provided to calculate the maximum impact of weathering on the carbon cycle. The temperature-sensitivity of the thermodynamically-limited silicate weathering provides a potential positive feedback to the carbon cycle which may shift the inner edge of the habitable zone.


INTRODUCTION
More than 20 Earth-sized planets 1 have been discovered in the habitable zones of their host stars. The classical habitable zone as used in studies such as Kasting et al. (1993) and Kopparapu et al. (2013) is the range of orbital distances, where an Earth-like planet can sustain liquid water on its surface for billions of years. The usual assumptions for the classical habitable zone include an atmosphere composed of CO 2 , H 2 O, and N 2 , as well 1 http://phl.upr.edu/projects/habitable-exoplanets-catalog as an Earth-like oceanic water reservoir. Within the next decade, the detected number of such Earth-sized planets is expected to increase due to space-based transit survey missions, such as PLATO (PLAnetary Transits and Oscillation of stars, Rauer et al. 2014). Nextgeneration telescopes including the James Webb Telescope (JWST ), the Extremely Large Telescope (ELT ) and the Thirty Meter Telescope (TMT ), or perhaps ESA's M-class mission ARIEL (Atmospheric Remotesensing Infrared Exoplanet Large-survey), might enable the atmospheric characterization of some of these planets.
Greenhouse gases such as CO 2 are essential in retaining heat in the atmosphere of Earth, however an excess of CO 2 along with water vapor may lead to a runaway greenhouse similar to Venus (Kasting et al. 1993;Kopparapu et al. 2013). To regulate the amount of CO 2 in the atmosphere, processes such as the weathering of rocks, degassing and regassing are essential (Ebelmen 1845;Urey 1952). One of the basic assumptions for the definition of the classical habitable zone is therefore the presence of the carbonate-silicate cycle (carbon cycle) that regulates the long-term climate (Kasting et al. 1993). However, it remains unclear how the carbon cycle may operate on rocky exoplanets where surface conditions could depart from those on modern Earth. The essence of the carbon cycle is captured by the Ebelman-Urey reaction involving the conversion of a silicate mineral (e.g., wollastonite CaSiO 3 ) to a carbonate mineral (e.g., calcite CaCO 3 ) in the presence of atmospheric CO 2 , CaSiO 3 + CO 2 ← → CaCO 3 + SiO 2 . (1) The reverse of reaction (1), metamorphism, converts carbonates back to silicates and releases CO 2 to the atmosphere. An important feature of the carbon cycle on Earth is the negative feedback of silicate weathering (e.g., Walker et al. 1981;Berner et al. 1983;Kump et al. 2000;Sleep & Zahnle 2001;Abbot et al. 2012;Foley 2015;Krissansen-Totton & Catling 2017;Graham & Pierrehumbert 2020). This feedback buffers the climate against changes in stellar luminosity and impacts the extent of the habitable zone (e.g., Kasting et al. 1993;Kopparapu et al. 2013). This feedback is sensitive to partial pressure of CO 2 in the atmosphere, surface temperature, and runoff (flow rate of water through soils) which facilitates weathering reactions and transport of aqueous chemical species from continents to oceans (e.g., Walker et al. 1981;Berner et al. 1983). High concentrations of CO 2 in the atmosphere elevate the surface temperature. High temperatures give rise to high evaporation and high precipitation rates. As precipitation intensifies, runoff also intensifies. As a result, silicate weathering intensifies, reducing the abundance of atmospheric CO 2 which in turn decreases the temperature. When temperatures become low, evaporation and precipitation rates become low thereby decreasing runoff and weathering. In the absence of intense weathering, volcanic degassing increases the abundance of CO 2 in the atmosphere sufficiently to increase the temperature again. In the absence of a self-regulating mechanism, the fate of Earth's atmosphere may have been like that of Venus (Gaillard & Scaillet 2014).
In addition to silicate weathering on continents, silicate weathering on the seafloor also has the potential to provide an equivalent negative feedback (e.g., Fran-cois & Walker 1992;Brady & Gíslason 1997;Sleep & Zahnle 2001;Coogan & Gillis 2013;Krissansen-Totton & Catling 2017;Charnay et al. 2017). Seafloor weathering is the low-temperature carbonation of the basaltrich oceanic crust facilitated by the circulation of seawater through hydrothermal systems. Seafloor weathering reactions, like continental weathering reactions, dissolve silicate minerals constituting rocks in the presence of water and CO 2 . Most studies modeling the carbon cycle expect the contribution of seafloor weathering to be smaller than continental weathering by up to a few orders of magnitude and thus neglect it (e.g., Walker et al. 1981;Berner et al. 1983;Caldeira 1995;Berner & Kothavala 2001). However, recent studies have shown that seafloor weathering is likely a significant component of global weathering during Earth's history (e.g., Coogan & Gillis 2013;Charnay et al. 2017;Krissansen-Totton et al. 2018) and is of critical importance on oceanworlds (e.g., Abbot et al. 2012;Foley 2015;Höning et al. 2019).
A prevalent assumption among studies modeling silicate weathering is that kinetics of mineral dissolution reactions determines the weathering flux where the subscript '0' denotes reference values, x CO2(g) is the volume mixing ratio of CO 2 in the gaseous state given by x CO2(g) = P CO2 /P with P CO2 as the CO 2 partial pressure and P as the total surface pressure, and β is constrained empirically from laboratory or field data (e.g., Walker et al. 1981;Berner et al. 1983;Kasting et al. 1993;Sleep & Zahnle 2001;Foley 2015;Krissansen-Totton & Catling 2017). The activation energy E a is also determined empirically from the dependence of kinetic rate coefficients on the Arrhenius law with T as the surface temperature and R as the universal gas constant (Palandri & Kharaka 2004;Brantley et al. 2008). Walker et al. (1981) adopt β = 0.3 based on kinetic rate measurements of feldspar weathering in a laboratory (Lagache 1965). Later studies adjusted the value of β using more laboratory and field measurements or balancing carbon fluxes (Brantley et al. 2008). A Bayesian inversion study performed by Krissansen-Totton & Catling (2017) for the carbon cycle on Earth based on data from past 100 Myr shows that the value of β for continental weathering is largely unconstrained: between 0.21−0.48 (prior: 0.2−0.5) for their nominal model and between 0.05−0.95 (prior: 0−1) for their Michaelis-Menten model (Volk 1987;Berner 2004). Moreover, the power-law exponent for seafloor weathering in Equation (2) is either assumed to be equal to 0.23 (Brady & Gíslason 1997) or varied between 0−1 (e.g., Sleep & Zahnle 2001;Coogan & Dosso 2015). More sophisticated formulations of seafloor weathering assume a dependence on the oceanic crustal production rate with yet another power-law exponent between 0−2 (Krissansen- Totton et al. 2018).
The flow of fluids such as rainwater through the porespace of soils facilitates silicate weathering reactions on continents. If the fluid spends less time in contact with the bedrock than the time needed for a reaction to attain chemical equilibrium, the reaction becomes rate limiting and reaction kinetics govern the amount of solute released (e.g., Aagaard & Helgeson 1982;Helgeson et al. 1984). The weathering is kinetically-limited. In contrast, at long fluid residence times, silicate weathering reactions reach chemical equilibrium resulting in a thermodynamic limit (also known as transport or runoff limits, Thompson 1959;Kump et al. 2000;Stallard & Edmond 1983;Maher 2011). For thermodynamicallylimited weathering, thermodynamic β values can be calculated from the reaction stoichiometry and equilibrium constants of weathering reactions. Winnick & Maher (2018) provide thermodynamic β values for the weathering of feldspar minerals (0.25−0.67). Observations that regional weathering fluxes on Earth depend strongly on runoff (e.g., Riebe et al. 2004) in addition to kinetics suggest that global weathering fluxes are a mixture of kinetically-and thermodynamically-limited regimes (Kump et al. 2000;Krissansen-Totton & Catling 2017). The transport-controlled approach (Maher 2010;Maher & Chamberlain 2014) allows the modeling of both thermodynamic and kinetic limits of weathering on temperate planets (Winnick & Maher 2018;Graham & Pierrehumbert 2020).
The lithology (rock type) is anticipated to play a role in the intensity of the silicate weathering feedback (e.g., Walker et al. 1981;Stallard & Edmond 1983;Kump et al. 2000). Winnick & Maher (2018) demonstrate that minerals in the feldspar mineral group exhibit weak to strong feedback at the thermodynamic limit of weathering. The kinetic weathering models are based on reaction rates of feldspars, common silicate minerals comprising granitic rocks on modern continents (e.g., Walker et al. 1981). However, rocks on the Hadean and Archean Earth were mafic (poorer in SiO 2 and feldspar content). The modern seafloor mostly contains basaltic rocks. Observations of clay minerals on Mars suggest past weathering processes (Bristow et al. 2018). Clues from stellar elemental abundances point to diverse planetary compositions (e.g., Wang et al. 2019;Spaargaren et al. 2020). Consequently, the surface lithologies (and therefore weathering fluxes) on rocky exoplanets are not necessarily similar to modern Earth.
In this study, we develop a silicate weathering model based on a transport-controlled approach (Maher & Chamberlain 2014) by considering fluid-rock reactions at chemical equilibrium for major silicate lithologies: CHILI (CHemical weatherIng model based on LIthology). This model tracks the total dissolved inorganic carbon in the aqueous surface reservoir assimilating four CO 2 drawdown regimes. In addition to continental weathering, this model is applied to seafloor weathering. The primary goal is to determine lithology-based weathering fluxes on the surface of temperate exoplanets by mitigating the impact of present-day Earth calibrations.

Proxies for Weathering
Continental weathering occurs in continental soils where runoff (flow rate of water through soils) facilitates fluid-rock weathering reactions ( Figure 1). Seafloor weathering occurs in the pore-space of the oceanic crust. Analogous to runoff, the fluid flow in off-axis lowtemperature hydrothermal systems facilitates seafloor weathering reactions (Stein & Stein 1994;Coogan & Gillis 2018). The fluid flow rate q (either runoff or hydrothermal fluid flow rate) is a free parameter in CHILI (Table 1). On Earth, the present-day continental runoff varies between 0.01−3 m yr −1 with a mean value of approximately 0.3 m yr −1 (Gaillardet et al. 1999;Fekete et al. 2002). The seafloor fluid flow rates have been measured between 0.001−0.7 m yr −1 with a mean value of about 0.05 m yr −1 (Stein & Stein 1994;Johnson & Pruis 2003;Hasterok 2013).
In the transport-controlled model of continental weathering, the weathering flux w is the product of the concentration of a solute of interest C and runoff q (e.g., Maher & Chamberlain 2014), We also apply the transport-controlled weathering model to seafloor weathering. The total silicate weathering rate 2 W is the sum of the continental (W cont ) and seafloor weathering rates (W seaf ). These rates are given by product of the continental (w cont ) or seafloor weathering flux (w seaf ) and the continental ( Figure 1. Silicate weathering model based on weathering reactions and fluid flow rates. The energy balance between incident stellar flux and outgoing longwave radiation results in evaporation and precipitation leading to runoff through the continental soils. Analogously, hydrothermal heat flux is responsible for the flow rate of hydrothermal circulation through the seafloor pore-space. In the transport-controlled model, the continental runoff and the seafloor fluid flow rate enable weathering reactions in the continental and seafloor pore-space, respectively. Key parameters and processes are labeled (see Tables 1, 3 and 2 for all parameters and output quantities). area and f is the continental area fraction: (4) To model silicate weathering, the aqueous bicarbonate ion concentration [HCO − 3 ] is normally used as a proxy for weathering (e.g., Winnick & Maher 2018). This is a good assumption at modern Earth conditions since HCO − 3 is the primary CO 2 -rich product of continental silicate weathering that is carried to oceans by rivers and reacts with Ca 2+ to precipitate Ca-rich carbonates on the seafloor. We are interested in modeling weathering at conditions more diverse than those on modern Earth. In highly alkaline conditions, the carbonate ion CO 2− 3 is produced in amounts similar to HCO − 3 . Moreover, aqueous carbon dioxide CO 2 (aq), although not produced by rock weathering, along with HCO − 3 and CO 2− 3 has the potential to contribute to Ca-Mg-Fe carbonate precipitation on the seafloor by making the solution acidic (e.g., calcite precipitation, Plummer et al. 1978). Together with HCO − 3 , CO 2− 3 and CO 2 (aq) give all of CO 2 drawdown from the atmosphere to the aqueous reservoir. Therefore, besides [HCO − 3 ] as a proxy for weathering (C = [HCO − 3 ] in Equation 3), we choose the total dissolved inorganic carbon (DIC) as a proxy for the total CO 2 drawdown (C = DIC in Equation 3) where

Major Silicate Lithologies
Previous studies modeling the transport-controlled continental weathering limit the lithology to individual minerals or one type of rock (e.g., feldspars or granite, Maher & Chamberlain 2014;Winnick & Maher 2018;Graham & Pierrehumbert 2020). Based on solar system observations, we consider three common silicate rocks that are exposed at the surface of a temperate planet: peridotite, basalt and granite. Peridotites are mantle rocks rich in MgO and FeO and poor in SiO 2 relative to basalts and granites. Basalts consist of minerals with Note-T is not a control parameter for some calculations (T = f (x CO 2 (g) ), Section 2.5).
relatively higher SiO 2 content. These igneous rocks are common examples of rocks present in the oceanic crust of Earth, lunar mare on the Moon, and the crust of Mars. Granites form the present-day continental crust on Earth. These are highly differentiated rocks that are rich in SiO 2 , Na 2 O and K 2 O due to partial melting and crystallization processes. Each of the three rocks are composed of 2−3 major mineral groups (Figure 2a). Each mineral group is similarly defined by 1−3 endmember minerals (Figure 2b). For instance, olivine is a solid solution of two endmember minerals, forsterite and fayalite. In the reduced set of mineral groups with only endmember minerals, we represent peridotite using pyroxenes (wollastonite and enstatite) and olivines (forsterite and fayalite), basalt using plagioclase feldspars (anorthite and albite) and pyroxenes (wollastonite, enstatite and ferrosilite), and granite using alkali feldspars (K-feldspar and albite), quartz and biotite micas (phlogopite and annite). The choice of endmember minerals makes these rocks idealized compared to natural rocks. Moreover, for a rock, not all endmember minerals of a mineral group can be modeled due to the unity activity assumption for endmember minerals (Section 2.3). Only major minerals in typical rock types are considered and minor minerals such as magnetite, hematite and pyrite are not discussed. Since the net contribution of carbonate weathering on the carbon cycle is small on timescales of the order of 100 kyr for Earth (e.g., Sleep & Zahnle 2001), the weathering of carbonate minerals is neglected.

Maximum Weathering Model
Concentrations of products of weathering reactions at chemical equilibrium allows to calculate a thermodynamic upper limit to weathering (maximum weathering). To calculate concentrations of silicate weathering products, a number of chemical reactions need to be considered. Atmospheric CO 2 (g) dissolves in water to produce CO 2 (aq), where g and aq represent gaseous and aqueous phases, respectively (Appendix C). Reactions between water, CO 2 (g) and silicate minerals produce the bicarbonate ion HCO − 3 . The bicarbonate ion further dissociates into the carbonate ion CO 2− 3 and H + . Water dissociates into H + and OH − . The relations between equilibrium constants of these reactions and thermodynamic activities of reactants and products are given in Appendix B. Thermodynamic activities are used to quantify the energetics of mixing of constituent components in solid or aqueous solutions (DeVoe 2001). Since no solid-solution behaviors are considered, the endmember minerals in this study involves no mixing, and their activities are thus unity. Furthermore, in a dilute solution like in this study, the H 2 O(aq) component is approximated by pure liquid water whose activity is also unity.
The activity of an aqueous species (e.g., CO 2 (aq)) is given by the product of its concentration [CO 2 (aq)] and the activity coefficient γ normalized to the standard state of concentration (DeVoe 2001). In an ideal, dilute solution, γ → 1 and [CO 2 (aq)] = a CO2(aq) C 0 , an assumption made throughout the study for all aqueous species. The activity of a gaseous species such as CO 2 (g) is given by the ratio of its fugacity in the gas mixture to the fugacity of pure CO 2 (g) at the total pressure P (i.e., surface pressure in this study). Thus, a CO2(g) = where P CO2 is the CO 2 partial pressure, P is the total surface pressure, Γ and Γ tot are fugacity coefficients for the CO 2 component and pure CO 2 , respectively. The fugacity coefficients give a correction factor for the nonideal behavior due to mixing and/or pressure effects. For CO 2 (g), the fugacity coefficient varies between 0.5 and 1 up to pressures of 200 bar (Spycher & Reed 1988). Our assumption of unity fugacity coefficients leads to a CO2(g) = PCO 2 P = x CO2(g) , where x CO2(g) is the CO 2 volume mixing ratio or the mole fraction of CO 2 (g). In this study, x CO2(g) is used as a control parameter instead of P CO2 since x CO2(g) allows us to probe weathering independently of the surface pressure.
Equilibrium constants depend on temperature and pressure (see Appendix A). Chemical reactions on continents are characterized by the surface temperature T and surface pressure P (Table 1). Seafloor weathering reactions are characterized by the seafloor pore-space temperature T and pressure P . In the temperature range 273−373 K, data suggest that T is within 1% of T (Krissansen-Totton & Catling 2017, and references therein) and hence T = T is assumed. Moreover, pressure in the range of 0.01−1000 bar has a negligible effect on equilibrium constants and resulting concentrations (see Figure A.1 and Figure B.3). Nonetheless, P is fixed to 200 bar which is the pressure at approximately 4 km depth in Earth's present-day oceans.
To demonstrate the computation of thermodynamic solute concentrations resulting from fluid-rock weathering reactions, weathering of peridotite is considered. Since two pyroxene endmembers (wollastonite and enstatite) and two olivine endmembers (forsterite and fayalite) are considered to constitute peridotite, dissolution reactions of these four minerals are of interest (Appendix B, Table B.1, rows: (a), (b), (d), (e)). Besides, reactions in the water-bicarbonate system are needed (Appendix B, Table B.1, rows: (o), (p), (q), (r)). These eight chemical reactions result in eight equations and ten unknowns. The eight equations are given by the relations between activities and equilibrium constants (Appendix B). The ten unknowns are the activities of CO 2 (g), CO 2 (aq), SiO 2 (aq), Ca 2+ , Mg 2+ , Fe 2+ , H + , OH − , HCO − 3 and CO 2− 3 . An additional equation is given by balancing the charges of all cations and anions present in the solution. These nine equations are further reduced to one polynomial equation with two unknowns, activities of HCO − 3 and CO 2 (g) (see Appendix B, Table B.2). If ferrosilite, the third endmember of pyroxene is also considered, an additional equation is introduced but the unknowns remain the same and the system is over-determined. This scenario stems from the assumption of endmember minerals with unity activities instead of solid solutions. In the case of solid solutions, the ac-tivities of endmember minerals become unknowns rather than being fixed at unity. This increases the number of unknowns for a given set of equations, and requires more equations that can be derived from additional reactions (e.g., see Galvez et al. 2015, for a treatment of solid solutions). Although mineral solid solutions are commonplace in rocks, endmember considerations allow us to establish a simple framework in which the effects of various lithologies can be explored on the weathering of exoplanets. Therefore, the presence of ferrosilite in peridotite is ignored. Moreover, the presence of Kfeldspar in basalt and anorthite in granite is ignored. For a given x CO2(g) , a HCO − 3 is calculated by finding the sole physical root of such a polynomial equation. The bicarbonate ion concentration at chemical equilibrium is obtained from the standard concentration, [HCO − 3 ] eq = a HCO − 3 × 1 mol dm −3 . Similarly, activities of all aqueous species are converted to concentrations using the standard concentration. Subsequently, all unknowns are calculated from the relations between activities and equilibrium constants (Table B.1).
At chemical equilibrium, the concentrations of DIC components increase and the pH of solution decreases with an increase in x CO2(g) at T = 288 K (Figure 3a). Since the activity of CO 2 (aq) is the product of K CO2 and x CO2(g) (reaction (o), Table B.1, [CO 2 (aq)] increases monotonically with x CO2(g) , where K CO2 as a function of T and P is obtained from the thermodynamic database (CHNOSZ, Dick 2019). The T -dependent relation between [CO 2 (aq)] and P CO2 provided by Pierrehumbert (2010, Equation 8.14) is accurate for a small range of T around 298 K (see Appendix C). In Figure 3(a), as x CO2(g) increases, [HCO − 3 ] eq and [CO 2− 3 ] eq increase monotonically, thereby increasing DIC eq , whereas pH eq decreases monotonically. For all except x CO2(g) < 20 ppmv, [HCO − 3 ] eq contributes significantly to DIC eq . The weathering flux (Equation 3) resulting from thermodynamic solute concentrations gives an upper limit to weathering (maximum weathering flux). Figure 3b shows that the thermodynamic weathering flux calculations using [HCO − 3 ] eq and DIC eq increase monotonically with x CO2(g) . The difference between the thermodynamic fluxes of [HCO − 3 ] eq and DIC eq at lower x CO2(g) arises from the excess contribution of [CO 2− 3 ] eq to DIC eq at conditions where the solution pH is basic. We follow the same procedure as peridotite to compute the maximum solute concentrations and the maximum weathering flux for basalt, granite or individual endmember silicate minerals by solving polynomial equations given in Appendix B and corresponding charge balance equations.

Generalized Weathering Model
In natural environments, fluid-rock reactions may not have the time to attain chemical equilibrium due to high fluid flow rates. The maximum weathering model (Section 2.3) does not accurately represent the weathering flux in the absence of chemical equilibrium for mineral dissolution reactions. We present the generalized weathering model by extending the transportcontrolled approach of Maher & Chamberlain (2014). Mineral dissolution reactions, being rate limiting, es-sentially buffer the concentration of reaction products including . In this limit, kinetics plays a more dominant role than thermodynamics. Maher & Chamberlain (2014) provide a solute transport equation to calculate such a transport-buffered (dilute) solute concentration from its value at chemical equilibrium, fluid flow rate (runoff) q and the Damköhler coefficient D w which gives the 'net reaction rate' (see below). The solute transport equation with HCO − 3 is In their formulation, Maher & Chamberlain (2014) multiply D w by an arbitrary scaling constant τ = e 2 ≈ 7.389 in order to scale up the solute concentration to 88% of its equilibrium value at q = D w . Ibarra et al. (2016) when applying the solute transport equation of Maher & Chamberlain (2014), use τ = 1 without any scaling. For simplicity, we use τ = 1, implying [HCO − 3 ] = 0.5 [HCO − 3 ] eq at q = D w . The resulting [HCO − 3 ] of our model at q = D w is about 43% lower than that of the Maher & Chamberlain (2014) model, which is a much smaller difference than the 5−10 orders of magnitude range of solute concentrations explored in this study.
The quantity Dw q in Equation (6) is the Damköhler number, the ratio of fluid residence time and chemical equilibrium timescale (Steefel & Maher 2009). The Dw q ratio is essentially the ratio of the 'net reaction rate' and the fluid flow rate. When the fluid residence time exceeds the chemical equilibrium timescale, or the net reaction rate exceeds the fluid flow rate (q < D w ), [HCO − 3 ] → [HCO − 3 ] eq and the weathering flux reaches its maximum value for a given q, w = [HCO − 3 ] eq q ( Figure 4). This weathering regime is the thermodynamically-limited weathering (hereafter, thermodynamic regime), also known as transport-or runoff-limited weathering. When the chemical equilibrium timescale exceeds the fluid residence time, or the fluid flow rate exceeds the net reaction rate (q > D w ), Dw q and w = [HCO − 3 ] eq D w , making w independent of [HCO − 3 ] eq because D w is modeled to be inversely proportional to [HCO − 3 ] eq (see below). This regime is known as kinetically-limited weathering (hereafter, kinetic regime). The transition between these two regimes occurs at q = D w . The comparison of timescales described here is conceptually identical to the 'quenching approximation' employed in atmospheric chemistry (Prinn & Barshay 1977;Tsai et al. 2017).
Rewriting the formulation of D w from Maher & Chamberlain (2014), where [HCO − 3 ] eq is the equilibrium solute concentration, k eff is the effective kinetic rate coefficient given by kinetics data (Appendix A), t s is the age of soils or porespace, A sp is the specific reactive surface area per unit mass of the rock, m is the mean molar mass of the rock, ψ = L(1 − φ)ρX r A sp is a dimensionless pore-space parameter that combines five parameters including flowpath length L, porosity φ, rock density ρ and fraction of reactive minerals in fresh rock X r . Although Maher & Chamberlain (2014) parameterize D w using nine quantities, in Appendix D, we show that D w is mostly sensitive to four of these quantities given their plausible range: [HCO − 3 ] eq , k eff , t s and L (L is absorbed in ψ). The remaining five parameters are fixed to reference values (Table 3). The parameter ψ scales k eff with dimensions of moles per unit reactive surface area of rocks per unit time to w with dimensions of moles per unit exposed continental or seafloor area for a given lithology per unit time. In Equation (7), the solid mass to fluid volume ratio ρ sf given in the D w formulation of Maher & Chamberlain (2014) is rewritten in terms of solid density and porosity using ρ sf = ρ(1 − φ)/φ. The kinetic rate coefficients of mineral dissolution reactions are obtained from Palandri & Kharaka (2004) as a function of T and pH (see Appendix A). The solution pH at chemical equilibrium is used to calculate k eff . For minerals with no k eff data, k eff of a corresponding endmember mineral from the same mineral group is adopted. For rocks, the minimum k eff among constituent minerals is used since the slowest reaction is rate limiting.
In Equation (7), when t s = 0 (young soils), D w = k eff ψ [HCO − ψ mAspts . This regime maybe termed as the 'slow kinetic' regime, however an accepted terminology is supply-limited weathering (hereafter, supply regime) which is limited by the supply of fresh rocks in the weathering zone (Riebe et al. 2003;West et al. 2005). An increase in t s decreases the influence of chemical kinetics on weathering. Although t s depends on the soil production and physical erosion rates which in turn are sensitive to climate, topography and fluid flow properties, there is no consensus on the formulation of soil production and physical erosion rates (Riebe et al. 2003;West et al. 2005;Gabet & Mudd 2009;West 2012;Maher & Chamberlain 2014;Foley 2015). Treating t s as a free parameter makes it possible to model both the age of continental soils and seafloor pore-space. The transition between the kinetic and supply regimes can be defined at t s = 1 k eff mAsp where the kinetic reaction rate equals the supply rate of fresh rocks (see Equation 7).
In the transport-controlled model, the mineral dissolution reactions (Appendix B, rows (a−n) in Table B.1) do not attain chemical equilibrium and their reaction products are given by the solute transport equation (Equation 6). On the other hand, reactions in the water-bicarbonate system (Appendix B, rows (a−n) in Table B.1), are considered to reach chemical equilibrium. This assumption is justified because the equilibrium timescale of chemical reactions in the waterbicarbonate system is of the order of milliseconds to days, unlike mineral dissolution reactions that may take months to thousands of years (e.g., Palandri & Kharaka 2004;Schulz et al. 2006). Thus, the generalized concentrations of aqueous species in the water-bicarbonate system such as [CO 3 ], x CO2(g) and equilibrium constants in Equation (5), DIC components and the solution pH for the generalized model of peridotite weathering as a function of x CO2(g) at a constant surface temperature (T = 288 K) are shown in Figure 4(a). Note that [CO 2 (aq)] is the same for the maximum and generalized cases (Figs. 3a 10 −1 10 1 10 3 10 5 Thermodynamic Kinetic CO 2 Diss.
(a) Generalized concentrations of DIC components and the solution pH of peridotite weathering as a function of x CO 2 (g) at T = 288 K and modern mean runoff of q = 0.3 m yr −1 (other parameters take reference values). From left to right, colored disks mark the transition between thermodynamic and kinetic regimes, and the inverted triangle marks the transition between kinetic and CO2 dissolution regimes with respect to DIC. (b) The corresponding flux of DIC and [HCO − 3 ]. The scales on vertical and horizontal axes are equal to those in Figure 3. and 4a) because [CO 2 (aq)] depends only on x CO2(g) and K CO2 (Table B.1). For x CO2(g) < 1 ppmv, all DIC components and pH in Figure 4(a) show the same behavior as in Figure 3(a) because q < D w , reiterating that chemical equilibrium calculations are valid in the thermodynamic regime resulting in the maximum weathering flux. The maximum ([HCO − 3 ] eq and DIC eq , Figure 3b) and generalized ([HCO − 3 ] and DIC, Figure 4b) weathering fluxes diverge from each other beyond the thermodynamic to kinetic regime transition that occurs at x CO2(g) = 1 ppmv. The choice of q determines this transition. At a higher q, the regime transition would shift to a lower x CO2(g) and vice versa. For x CO2(g) > 1 ppmv, [HCO − 3 ] becomes independent of x CO2(g) in the kinetic weathering regime. This is because the kinetic rate coefficient of the slowest reaction (fayalite dissolution) is constant at a fixed T and does not vary when pH is basic (see Appendix A). However, the results shown in Figure 4 are expected to change when temperature is not held constant and depends on x CO2(g) where the relation between the two is given by a climate model (Section 2.5).
In the kinetic weathering regime, shows a strong decrease because of its sole dependence on x CO2(g) . For x CO2(g) > 23 000 ppmv, [CO 2 (aq)] contributes significantly to DIC ( Figure 3a) as a result of the CO 2 dissolution reaction (reaction (o) in Table B.1). Although [CO 2 (aq)] is not a product of rock weathering, [CO 2 (aq)] has the potential to contribute to carbonate precipitation on the seafloor by producing acidic conditions (see Section 2.1). Together with the three weathering regimes (thermodynamic, kinetic and supply), the CO 2 dissolution regime constitutes the fourth CO 2 drawdown regime. We define the transition between the kinetic to CO 2 dissolution regimes at , which is half of the total DIC. We follow the same procedure as peridotite to compute the generalized concentrations for the weathering of basalt, granite or individual minerals and find that lithology strongly impacts the occurrence of CO 2 drawdown regimes (Appendix B, Fig. B.2).

Climate Model
The greenhouse effect of CO 2 exerts a strong control on the planetary surface temperature. A climate model enables us to express the surface temperature T as a function of the CO 2 volume mixing ratio x CO2(g) = P CO2 /P , where P CO2 is the CO 2 partial pressure and P is the surface pressure for a given planetary albedo α and top-of-atmosphere stellar flux S. Such a climate model is essential to assess the role of climate in weathering on temperate planets. CHILI provides the functionality to couple T and x CO2(g) using any climate model. Previous studies provide formulations of climate models (e.g., Walker et al. 1981;Kasting et al. 1993). Recently, Kopparapu et al. (2013Kopparapu et al. ( , 2014 performed 1D radiative-convective calculations to obtain T as a function of P CO2 , α and S. Studies such as Haqq-Misra et al. (2016); Kadoya & Tajika (2019) provide fitting functions to the models of Kopparapu et al. (2013Kopparapu et al. ( , 2014. We use the fitting function provided by Kadoya & Tajika (2019) to couple T in the range 150−350 K with P CO2 in the range 10 −5 − 10 bar at saturated H 2 O and 1 bar N 2 . For α = 0.3 (present-day albedo of Earth) and S = 1360 W m −2 (present-day solar flux), the Kadoya & Tajika (2019) fitting function results in T between 280−350 K for x CO2(g) between 10 ppmv and 2 × 10 5 ppmv (see Appendix E for details).

Maximum Weathering for Various Lithologies
The thermodynamic solute concentrations provide an upper limit to weathering (e.g., peridotite weathering, Figure 3). To evaluate the case of maximum weathering on temperate planets, in Figure 5, we provide the [HCO − 3 ] eq weathering flux of three common rocks (basalt, peridotite, granite) as a function of climate properties, x CO2(g) and T , at present-day mean runoff of 0.3 m yr −1 . The dependence of weathering on total surface pressure is negligible (see Appendix B). The [HCO − 3 ] eq weathering flux for all rocks increases monotonically with x CO2(g) at a constant temperature because weathering intensifies as x CO2(g) increases (Figure 5a). This is a direct consequence of the calculations of solute concentrations at chemical equilibrium (Section 2.3). Table 4 gives fitting parameters of the thermodynamic weathering flux of [HCO − 3 ] eq to the kinetic weathering expression (Equation 2) for silicate rocks and minerals considered in this study. Such a fit, although not the best approximation of the calculated values, provides a way to compare the sensitivity of thermodynamic weathering to climate properties with studies assuming kinetic weathering (e.g., Walker et al. 1981;Sleep & Zahnle 2001;Foley 2015).
Fluid-rock reactions produce both monovalent (e.g., K + , Na + ) and divalent cations (e.g., Ca 2+ , Mg 2+ ). The sensitivity of [HCO − 3 ] eq flux to x CO2(g) depends on the capacity of a rock to produce divalent cations. More divalent cations in the solution require more HCO − 3 to balance the charges. Consequently, the higher the fraction of divalent cations in a rock, the higher the thermodynamic x CO2(g) sensitivity (β). For instance, peridotite, which produces only divalent cations, exhibits the highest β among the three rocks ( Figure 3). Granite produces more monovalent cations than peridotite and basalt, and therefore has the lowest β. This effect 10 0 10 1 10 2 10 3 10 4 10 5 10 6  has been discussed by Ibarra et al. (2016); Winnick & Maher (2018). Even with a smaller β, basalt shows a higher thermodynamic weathering flux than peridotite because of the presence of highly-weatherable anorthite.
The β values for endmember minerals within the same mineral group (pyroxene, olivine, mica and amphibole), Note-w0 is not a fitting parameter but the value of w at x CO 2 (g) = 280 ppmv, T = 288 K and q = 0.3 m yr −1 .
with the exception of feldspars, are similar to each other (Table 4). This result is again attributed to the presence of monovalent or divalent cations. The divalent cation-producing pyroxene, olivine, mica and amphibole endmembers exhibit β ∼ 0.5 and monovalent cationproducing albite and K-feldspar exhibit β ∼ 0.25. The deviation from these ideal values of 0.5 and 0.25 is a result of the simultaneous consideration of the mineral dissolution reaction and the water-bicarbonate reactions. For instance, Winnick & Maher (2018) show that reaction stoichiometry controls β values by considering dissolution reactions of individual feldspar minerals and find β = 0.25 for albite and K-feldspar. However, we find that the presence of ions produced by the waterbicarbonate system makes β dependent on equilibrium constants of all reactions considered in addition to stoichiometry and hence β = 0.26 for albite and K-feldspar (Table 4). The β value of the divalent cation-producing anorthite is 0.72, considerably higher than other divalent cation-producing minerals. This is again attributed to reaction stoichiometry as demonstrated by Winnick & Maher (2018), although their study results in β = 0.67 because of neglecting reactions in the water-bicarbonate system. Contrary to observations of the increase in kinetic weathering flux with temperature (e.g., Walker et al. 1981), the thermodynamic weathering flux for rocks decreases with temperature at a fixed x CO2(g) (Figure 5b). The fitting parameter, the activation energy of thermodynamic weathering provides a scaling relation between weathering and temperature. Unlike kinetic weathering, thermodynamic E a is negative for weathering of all rocks and minerals except albite and K-feldspar (Table 4). This result is a consequence of the decrease in equilibrium constants as a function of temperature for all mineral dissolution reactions except those of albite and K-feldspar (see Appendix A). Winnick & Maher (2018) observe this effect for plagioclase feldspars which contain anorthite in addition to albite. This effect is discussed further in Section 4.2. The fitting parameters provided in Table 4 should be used with caution as β depends on T and E a depends on x CO2(g) . For example, β for peridotite at x CO2(g) = 280 ppmv varies between 0.71 at T = 273 K and 0.77 at T = 373 K. Whereas, E a for peridotite at T = 288 K varies between −55 kJ mol −1 at x CO2(g) = 1 ppmv and −46 kJ mol −1 at x CO2(g) = 10 6 ppmv. This provides further reason to calculate thermodynamic concentrations consistently by considering all necessary reactions simultaneously, as formulated in this study. Figure 6 shows the sensitivity of both maximum and generalized weathering fluxes of peridotite to CO 2 volume mixing ratio and surface temperature. The maximum weathering model gives an upper limit to the weathering flux. The generalized weathering flux is either equal to or smaller than the maximum weathering flux depending on if the generalized model encounters the thermodynamic regime or not. The generalized weathering flux in the kinetic regime is lower than that in the thermodynamic regime and higher than that in the supply regime. Since the weathering flux in the thermodynamic regime is an upper limit to weathering, the kinetic flux cannot exceed the thermodynamic flux. The CO 2 dissolution flux exceeds the weathering fluxes in supply and kinetic regimes at high CO 2 volume mixing ratios. To isolate the effect of x CO2(g) on weathering, we present peridotite weathering flux as a function of x CO2(g) at a fixed temperature (Figure 6a). The sen-sitivity of weathering to temperature is demonstrated at two x CO2(g) values (Figure 6b,d). In reality, T depends on x CO2(g) due to the greenhouse effect of CO 2 which is normally modeled using a climate model (Section 2.5). The strength of the coupling between T and x CO2(g) is sensitive to numerous parameters including incident solar flux and planetary albedo that are fixed to present-day values (Table 3). Since the solar flux is held constant, the coupling between T and x CO2(g) becomes stronger than that during Earth's history where the solar flux dropped to about 70% of its present-day value. Figure 6(c) demonstrates peridotite weathering under a limit of strong coupling between T and x CO2(g) .

Climate Sensitivity of Peridotite Weathering
The maximum (thermodynamic) DIC and [HCO − 3 ] fluxes increase monotonically with x CO2(g) at T = 288 K (Figure 6a). This monotonic behavior is a direct result of chemical equilibrium calculations with x CO2(g) as a free parameter at a fixed temperature. The difference between the maximum DIC and [HCO − 3 ] fluxes at low x CO2(g) is due to the excess contribution of [CO 2− 3 ] to DIC when the solution pH is basic (see Figure 3). The thermodynamic weathering flux decreases with T at a fixed x CO2(g) (Figure 6b,d). As explained in Section 3.1 and later discussed in Section 4.2, this decrease is a result of the decrease in equilibrium constants of mineral dissolution reactions as a function of temperature (see Appendix A). The thermodynamic flux at x CO2(g) = 10 5 ppmv is higher than at x CO2(g) = 280 ppmv by about a factor of 30. Unlike at x CO2(g) = 280 ppmv, at x CO2(g) = 10 5 ppmv there is a negligible difference between the thermodynamic DIC and [HCO − 3 ] fluxes. When surface temperature depends on x CO2(g) via a climate model, the behavior of thermodynamic weathering as a function of x CO2(g) depends on the trade off between the individual effects of T and x CO2(g) (Figure 6c). Up to x CO2(g) = 10 4 ppmv, the thermodynamic fluxes increase with x CO2(g) because the x CO2(g) effect dominates. Whereas, beyond x CO2(g) = 10 4 ppmv, there is a small decrease in the thermodynamic fluxes because the effect of T takes over.
function of T at constant x CO2(g) , kinetic weathering exhibits a strong dependence on temperature as seen in Figure 6(b,d). For the x CO2(g) = 280 ppmv case, the generalized (t s = 0) model switches from the ki-netic to thermodynamic regime at ∼310 K as the upper limit of weathering is encountered. For the x CO2(g) = 10 5 ppmv case, this transition temperature increases to ∼340 K. The transition from kinetic to thermody-namic regime occurs at high T and low x CO2(g) . The sequence of this transition is in contrast to the switch from the thermodynamic to kinetic regime that occurs at x CO2(g) = 1 ppmv as a function of x CO2(g) (see Figure 4). If the climate model is invoked, the generalized (t s = 0) weathering flux increases steeply with x CO2(g) and encounters the thermodynamic regime at about x CO2(g) = 10 5 ppmv, a result of the kinetic temperature dependence (Figure 6c). When T is strongly coupled to x CO2(g) , the generalized model may switch from the thermodynamic to kinetic regime at low x CO2(g) and from the kinetic regime back to the thermodynamic regime at high x CO2(g) . When soils are old (present-day mean soil age, t s = 100 kyr), the [HCO − 3 ] flux at T = 288 K as a function of x CO2(g) is constant and in the supply regime (Figure 6a). This regime is limited by the supply of fresh rocks and the influence of chemical kinetics on weathering is small compared to when soils are young (t s is small). This regime is almost independent of T (Figure 6b,d). Consequently, unlike the kinetic weathering flux, the supply-limited weathering flux is independent of x CO2(g) when T -dependence via the climate model is invoked (Figure 6c). In the supply regime, the weathering flux depends on the age of soils which is held constant. There is almost no difference between the generalized (t s = 100 kyr) DIC and [HCO − 3 ] fluxes between Figure 6(a) and (c) because the inclusion of the T effect via the climate model is not strong enough to escape the transition from the supply to CO 2 dissolution regime. In contrast, for the generalized (t s = 0) model, invoking the climate model increases the weathering fluxes in the kinetic regime enough to enter the thermodynamic regime at about x CO2(g) = 10 5 ppmv.
The generalized DIC fluxes for t s = 0 and t s = 100 kyr models enter the CO 2 dissolution regime at about 5000 ppmv and 2.5 × 10 4 ppmv, respectively (Figure 6a). This regime transition is not a result of rock weathering but rather due to excess atmospheric CO 2 dissolution in water bodies such that [CO 2 (aq)] > DIC/2 (Section 2.4). Since the [CO 2 (aq)] flux depends linearly on x CO2(g) , the slope of the DIC flux in the CO 2 dissolution regime is unity (Figure 6a,c). At low x CO2(g) , none of the generalized weathering models enter this regime (Figure 6b). At high x CO2(g) , the DIC flux at t s = 100 kyr is in the CO 2 dissolution regime for the entire T range and the DIC flux at t s = 0 is in the CO 2 dissolution regime for T < 297 K (Figure 6d). The CO 2 dissolution DIC flux decreases with T because the equilibrium constant of the CO 2 dissolution reaction also decreases with temperature (Appendix A).

Endmember Cases of Continental and Seafloor Weathering
We apply the generalized weathering model to both continental (Figure 7a,b) and seafloor silicate weathering (Figure 7c,d) for diverse cases that may represent weathering scenarios on temperate rocky exoplanets. The climate model is used to couple T with x CO2(g) by fixing the stellar flux and planetary albedo to presentday Earth values (Section 2.5). This implies a strong coupling between T with x CO2(g) such that T varies from 280 K to 350 K when x CO2(g) varies from 10 ppmv to 2 × 10 5 ppmv. The lithology and pore-space properties of continents and seafloor on present-day Earth are different from each other. This presents an opportunity to test the generalized weathering model in an extended parameter-space beyond applications to continental weathering. As basalt and granite represent the endmember cases of rock lithologies considered in this study, we show only the results for basalt and granite in Figure 7. The sensitivity of the DIC weathering flux to x CO2(g) is a complex function of climate, fluid flow rate, rock and pore-space properties. This result implies that the weathering flux cannot be simply approximated by Equation (2) assuming either kinetic weathering (e.g., Walker et al. 1981;Berner et al. 1983;Sleep & Zahnle 2001;Foley 2015) or thermodynamic weathering (Winnick & Maher 2018, and Section 3.1, this study).
In Figure 7(a), the continental DIC flux is calculated for three values of runoff that are representative of arid, modern-mean and humid climates. Gaillardet et al. (1999) report regional variations in the presentday runoff from 0.01−3 m yr −1 . Since the choice of humid runoff in our model is higher than the arid runoff by three orders of magnitude, the corresponding weathering fluxes should differ by the same amount if the weathering fluxes are in the thermodynamic regime. However, the differences are smaller in the given x CO2(g) range. This is because non-thermodynamic regimes exhibit smaller weathering fluxes than the thermodynamic regime. For example, up to x CO2(g) = 200 ppmv, the arid model of granite is in the thermodynamic regime (upper limit for this model), whereas the humid model of granite is in the supply regime (lower limit for this model). Thus, the differences between the two models are smaller than the maximum possible differences. For the arid model in Figure 7(a), basalt exhibits supply-limited weathering for the given x CO2(g) range, and granite is in the thermodynamic regime up to x CO2(g) = 2000 ppmv. The excess DIC flux for basalt at x CO2(g) < 100 ppmv arises due to the contribution of [CO 2− 3 ] to DIC. For more arid climates, lithology plays an even more important role as the weathering becomes transport-or 10 1 10 2 10 3 10 4 10 5 x CO 2 (g) [ppmv] , where T = f (x CO 2 (g) ) is given by the climate model (Section 2.5). Colored disks denote the transition between thermodynamic and kinetic/supply regimes. Inverted triangles denote the transition between kinetic/supply and CO2 dissolution regimes.
thermodynamically-limited for the whole x CO2(g) range.
In contrast, the supply and CO 2 dissolution regimes are independent of lithology. For this reason, lithology has no impact on the modern-mean and humid cases. The age of soils is a key parameter that determines if the weathering is limited by reaction kinetics or limited by the supply of fresh rocks. Figure 7(b) shows that lithology has no influence on the weathering flux of old soils (t s = 100 kyr and t s = 1 Myr) as opposed to young soils (t s = 0). The old soil models are either in the supply or CO 2 dissolution regimes. Granite and basalt in young soils are in different weathering regimes. Granite is in the thermodynamic regime for the whole x CO2(g) range because the net reaction rate (D w = 0.33−1 m yr −1 ) is higher than the fluid flow rate (q = 0.3 m yr −1 ). In contrast, D w (0.001−0.5 m yr −1 ) of basalt is lower than q = 0.3 m yr −1 for most of the x CO2(g) range, implying a transition from kinetic to thermodynamic above x CO2(g) ∼ 10 5 ppmv.
On Earth, the continental and seafloor pore-space differ in terms of pore-space properties besides the differences in lithology. The dimensionless pore-space parameter ψ depends on the flowpath length L that is normally assumed to be of the order of the regolith thickness (Section 2.4). Since the thickness of the oceanic crust where seafloor weathering occurs is of the order of 100 m (e.g., Alt et al. 1986;Coogan & Gillis 2018), L = 100 m and ψ = 100ψ 0 . Moreover, the average age of the oceanic crust on Earth at present day is approximately equal to 50 Myr, about 500 times the mean age of continental soils.
In Figure 7(d), we compare the seafloor DIC flux for modern-mean values of t s = 50 Myr and ψ = 100ψ 0 with two models, one with t s = 0 and another with ψ = ψ 0 . The modern-mean seafloor DIC flux is either in the supply or CO 2 dissolution regime making it independent of lithology. The supply-limited fluxes of seafloor weathering are smaller than those of continental weathering by a factor of 5 because w ∝ ψ ts in this regime. When ψ is lowered from 100ψ 0 to ψ 0 , the supply-limited weathering flux decreases by two orders of magnitude and the CO 2 dissolution limit is attained at a much lower x CO2(g) . For the young seafloor pore-space case, the DIC fluxes of basalt and granite are in the thermodynamic regime, and consequently the impact of lithology is pronounced. Compared to the t s = 0 basalt model in Figure 7(b) that is in the kinetic regime, the t s = 0 basalt model in Figure 7(d) is in the thermodynamic regime. This difference arises due to the choice of ψ that makes the net reaction rate (D w ) higher than the fluid flow rate, pushing basalt into the thermodynamic regime. Being in the thermodynamic regime, chemical equilibrium controls the weathering flux of basalt. The increase of the weathering flux of basalt with x CO2(g) at low x CO2(g) is not evident because of the decreasing contribution of [CO 2− 3 ] to DIC with increasing x CO2(g) . At high x CO2(g) , the thermodynamic weathering flux of basalt decreases as the decreasing effect of T takes over the increasing effect of x CO2(g) . This decreasing effect of T is due to the decrease in equilibrium constants of weathering reactions as a function of temperature (Appendix A). The effect of T on the weathering flux of granite is small and hence there is no net decrease in the weathering flux at high x CO2(g) .
In Figure 7(c), we calculate the seafloor DIC flux for three hydrothermal fluid flow rates, where the two extreme q values differ by three orders of magnitude, similar to the strategy in Figure 7(a). Since the fluid flow rates are directly proportional to the hydrothermal heat flux (Stein & Stein 1994;Coogan & Gillis 2013), a variation in the hydrothermal heat flux implies a variation in the fluid flow rate. Depending on the age of oceanic crust, present-day fluid flow rates are observed from 0.001−0.7 m yr −1 (Johnson & Pruis 2003). The three cases of seafloor weathering shown in Figure 7(c) are broadly similar to their continental counterparts. The modern-mean and high fluid flow rates result in lithology-independent weathering fluxes that are either in the supply-limited or CO 2 dissolution regimes. The low hydrothermal fluid flow rate causes the granite model to be in the thermodynamic regime up to x CO2(g) = 2 × 10 4 ppmv. The excess DIC flux for basalt at x CO2(g) < 500 ppmv arises due to the contribution of [CO 2− 3 ] to DIC.

Weathering Regimes and the Role of Lithology
In this study, concentrations of silicate weathering products are calculated by the simultaneous consideration of dissolution reactions of all minerals present in a rock as well as reactions in the water-bicarbonate system. Three common silicate rocks (peridotite, basalt and granite) are examined. We develop the maximum weathering model (Section 2.3) presuming chemical equilibrium and the generalized weathering model (Section 2.4) that applies to both equilibrium and nonequilibrium conditions. The generalized weathering model allows us to explore weathering in four CO 2 drawdown regimes (Figure 8). To simulate the transportbased dilution of equilibrium concentrations of weathering products, the solute transport equation of Maher & Chamberlain (2014) is implemented. This equation is based on the interplay between the fluid flow rate (q) and the net reaction rate (D w ). When q < D w , reactions are limited by thermodynamics or transport (runoff) and the weathering is thermodynamically-limited (Figure 8a). In contrast, when q > D w , reactions are limited by kinetics (and independent of runoff) and the weathering is kinetically-limited.
Another parameter, the age of soils (t s ), is introduced by Maher & Chamberlain (2014) to model the effect of limited supply of fresh rocks on the net reaction rate. The higher the t s , the lower is the D w . This divides the kinetic regime into another regime, supply-limited weathering (Figure 8a). Either at high runoff or high CO 2 volume mixing ratios, the contribution of [CO 2 (aq)] to DIC exceeds that of [HCO − 3 ] + [CO 2− 3 ]. At such conditions, a fourth CO 2 drawdown regime takes over: CO2 drawdown regimes with respect to DIC flux of rocks for two generalized weathering models (ts = 0 and ts = 100 kyr). (a) DIC flux as a function of fluid flow rate at x CO 2 (g) = 280 ppmv and T = 288 K. The vertical black line is the modern-mean runoff (q = 0.3 m yr −1 ). (b) DIC flux as a function of x CO 2 (g) at modern-mean runoff of q = 0.3 m yr −1 where T = f (x CO 2 (g) ) is given by the climate model (Section 2.5). The vertical black line is the pre-industrial CO2 volume mixing ratio (x CO 2 (g) = 280 ppmv). Colored disks mark the transition between thermodynamic and kinetic/supply regimes and colored inverted triangles denote the transition between kinetic/supply-limited and CO2 dissolution regimes.
CO 2 dissolution. Like the thermodynamic regime, the CO 2 dissolution regime is transport-limited (Figure 8a). Although [CO 2 (aq)] is not produced by rock weathering, together with [HCO − 3 ] and [CO 2− 3 ], [CO 2 (aq)] has the potential to contribute to carbonate precipitation on the seafloor in acidic conditions. Therefore, it is essential to track total dissolved inorganic carbon in the surficial aqueous reservoir.
As a function of transport and climate properties, the DIC flux shows a strong dependence on lithology in the thermodynamic and kinetic regimes and a weak dependence on lithology in the supply and CO 2 dissolution regimes (Figure 8). This impact of lithology for the four CO 2 drawdown regimes as a function of q is seen in Figure 8(a). There is about an order of magnitude difference between the thermodynamic fluxes of basalt and peridotite as well as those of peridotite and granite. The thermodynamic weathering flux is proportional to the equilibrium DIC concentration which is strongly sensitive to the mineralogy considered for a given rock. In the kinetic regime, the differences are smaller but significant. The kinetic weathering flux is proportional to the effective rate coefficient of a given rock. In contrast, in the supply-limited and CO 2 dissolution regimes, the DIC fluxes almost overlap with each other because lithology does not impact these regimes.
The supply regime is strongly sensitive to pore-space space properties that are held constant. The CO 2 dissolution regime depends only on climate and transport properties (w = K CO2 x CO2(g) q). The presence of CO 2 dissolution regime indicates that rock weathering is a less effective pathway to draw down CO 2 from the atmosphere to the surface.
The climate sensitivity of the weathering flux showcases that the weathering regime may depend on lithology. Figure 8(b) shows that the generalized t s = 0 model of the weathering fluxes of peridotite and basalt increase steeply with x CO2(g) up to ∼ 10 5 ppmv where T = f (x CO2(g) ) is obtained from the climate model. Beyond this x CO2(g) value, peridotite and basalt enter the thermodynamic regime. This is because the weathering flux in the kinetic regime cannot keep increasing indefinitely. As soon as the model hits the thermodynamic upper limit, the model follows the thermodynamic sensitivities of weathering that are independent of porespace properties and depend on chemical equilibrium and lithology. Important to note that extrapolations of kinetic weathering expressions (Equation 2) used in previous studies (e.g., Walker et al. 1981;Berner et al. 1983;Sleep & Zahnle 2001;Foley 2015; Krissansen-Totton & Catling 2017) may incorrectly predict the weathering flux to be higher than the upper limit provided by the thermodynamic flux.
Unlike basalt and peridotite, the generalized t s = 0 model of granite is in the thermodynamic regime for the given x CO2(g) range (Figure 8b). This is also evident from Figure 8(a) where the vertical line (q = 0.3 m yr −1 ) falls right before the thermodynamic to kinetic regime transition for granite. Thus, the climate sensitivity of granite weathering at q = 0.3 m yr −1 is determined largely by thermodynamics instead of kinetics.
On the other hand, for the generalized models at the present-day mean soil age of t s = 100 kyr, the weathering fluxes of the three rocks overlap with each other. This is because this t s value is so high that the effect of k eff on the 'net reaction rate' D w is negligible in pushing the model out of the supply regime. In this regime, the weathering flux is independent of the climate and transport properties. Beyond x CO2(g) ∼ 10 4 ppmv, the models are in the CO 2 dissolution regime, implying that lithology no longer dictates weathering.
Since most laboratory measurements of kinetic rate coefficients are available for individual minerals, previous studies discuss weathering of individual minerals instead of rocks (e.g., Walker et al. 1981;Berner et al. 1983). In reality, all minerals in rocks undergo weathering contemporaneously, rendering consideration of individual minerals in isolated systems as less informative. Since all minerals in a rock are in contact with the aqueous solution, solute concentrations are buffered by the dissolution reactions of these minerals. It is essential to consider these reactions simultaneously to solve for solute concentrations. The generalized weathering model shows that the choice of individual minerals or rocks determines the weathering regime. For example, the feldspar endmember minerals (anorthite, albite and K-feldpsar) are mostly in the thermodynamic regime for x CO2(g) between 100 ppmv and 10 4 ppmv, whereas rocks are mostly in the kinetic regime (see Figure B.2). In their implementation of the transport-controlled model, Graham & Pierrehumbert (2020) find that weathering is largely independent of kinetics because of their choice of oligoclase (a type of plagioclase feldspar mineral) to model weathering, which is in the thermodynamic regime of weathering for a wide range of CO 2 partial pressures similar to feldspar endmembers shown in Figure B.2. This study shows that inferring the weathering flux of a rock from only one of its constituent minerals is misleading.

Positive Feedback of Weathering at High Temperature
A common understanding of present-day weathering on Earth is that weathering intensifies with surface temperature (e.g., Walker et al. 1981;Berner et al. 1983;Kump et al. 2000;Brantley et al. 2008). This Tdependence of weathering is due to the increase in kinetic rate coefficients of mineral dissolution reactions as a function of T (see Appendix A). Laboratory and field measurements of kinetic rate coefficients are fitted to the Arrhenius law, w ∝ exp (−E a /RT ), where E a is the activation energy and R is the universal gas constant (Palandri & Kharaka 2004). The generalized weathering model based on the transport-controlled approach (Maher & Chamberlain 2014) captures the diversity of weathering regimes in a generic formulation, which is particularly useful for application to exoplanets with potentially diverse surface environments. Studies applying the transport-controlled weathering model find that T has a small effect on weathering (Winnick & Maher 2018;Graham & Pierrehumbert 2020). This statement holds in the thermodynamic regime of weathering for certain plagioclase feldspars which coincidentally exhibit a small T sensitivity, although Winnick & Maher (2018) find that the equilibrium [HCO − 3 ] resulting from plagioclase feldspars decreases with T .
To isolate the effect of temperature on weathering from that of the x CO2(g) effect, we first present weathering as a function of T at a fixed x CO2(g) in Figure 9(a). The generalized weathering model at t s = 0, q = 0.3 m yr −1 and x CO2(g) = 280 ppmv (no climate model) shows that the DIC flux of rocks is in the kinetic regime at low temperatures and in the thermodynamic regime at high temperatures. In the kinetic regime, there is a steep increase in the weathering flux of granite up to ∼285 K, peridotite up to ∼315 K, and basalt up to ∼325 K because the effective kinetic rate coefficients show an exponential increase with temperature. Beyond these transition temperatures (where the net reaction rate equals the fluid flow rate), the models enter their respective thermodynamic regimes. This transition implies a switch from the negative feedback of silicate weathering to the carbon cycle to the positive feedback. Since the thermodynamic regime gives the maximum possible weathering flux, the kinetic weathering flux cannot exceed the thermodynamically-limited weathering flux in the generalized model. Importantly, the kinetic weathering flux of granite is higher than that of basalt and peridotite because the slowest kinetic reaction among the constituent endmember minerals of granite has a higher kinetic rate coefficient than those of basalt and peridotite. This result is based on the choice of endmember minerals to define rocks and is ex- pected to change with the choice of other endmember minerals or mineral solid solutions. The kinetic weathering flux increases with T as expected, however the thermodynamic weathering flux decreases with T ( Figure 9a). As described in Section 3.1, if fitted to a kinetic weathering expression (2), the thermodynamic model gives a negative value for the activa-tion energy (E a ) for the weathering of rocks and most minerals. Such a negative value highlights that thermodynamic weathering exhibits a negative slope as a function of T . This negative exponential decrease in weathering is a result of the negative slope exhibited by equilibrium constants of mineral dissolution reactions as a function of T except for albite and K-feldspar (see Appendix A). This T -dependence is a thermodynamic property of mineral dissolution reactions in the aqueous system. For rocks, constituent minerals determine the overall dependence of thermodynamic weathering on T . For example, for a rock consisting of only albite and Kfeldpsar, the thermodynamic weathering flux is expected to increase with T , similar to the prevalent understanding of the T -dependence of kinetic weathering. However, apart from these two minerals, all other minerals considered in this study show a negative slope. Since granite consists of both albite and K-feldspar in addition to quartz, phlogopite and annite, granite shows the least negative slope among the three rocks considered (Figure 9a). Moreover, weathering of basalt in the thermodynamic regime decreases with T more steeply than peridotite because of the presence of anorthite in basalt that exhibits the most negative slope in the equilibrium constant−temperature parameter-space. When the climate model is invoked assuming presentday solar flux and present-day planetary albedo, x CO2(g) increases from 10 ppmv to 2 × 10 5 ppmv as a function of the T varying from 280 K to 350 K and the decreasing effect of T on weathering disappears (Figure 9b). This is because in the thermodynamic regime the weathering flux shows a positive power-law dependence on x CO2(g) , which is stronger than the exponential decrease in weathering with temperature. And in the kinetic regime, the exponential increase in weathering with temperature dominates. Interesting to note that basalt and peridotite enter the thermodynamic regime only at ∼340 K and ∼345 K respectively, whereas granite is in the thermodynamic regime for the given temperature range. For this climate model assuming a constant stellar flux, results in a strong coupling between T and x CO2(g) , providing another extreme as opposed to when x CO2(g) is held constant. Depending on the stellar flux and planetary albedo, the strength of coupling between T and x CO2(g) is probably in between the two panels shown in Figure 9, pushing the thermodynamic to kinetic transition temperature given in Figure 9(a) to higher values. The supply-limited weathering flux, on the other hand, is independent of T (e.g., Figure 6c,d). This is because the supply-limited weathering flux depends on pore-space parameters including t s and ψ that are assumed to be constant for a given model. In reality, the age of soils may indirectly depend on temperature through the effect of precipitation and runoff on physical erosion rates (e.g., West et al. 2005;West 2012;Foley 2015).
Previous studies have argued that since runoff depends on precipitation which in turn depends on T , the weathering flux should increase even more strongly with T when a T -dependent runoff is assumed (e.g., Berner & Kothavala 2001). Figure 9 shows that this statement does not hold in the kinetic regime of weathering. This is because the kinetic regime is independent of runoff in the transport-controlled model (Maher & Chamberlain 2014). Considering a linear dependence of runoff on temperature (e.g., Equations 41, 42, Graham & Pierrehumbert 2020) instead of a constant runoff, there is no impact on weathering because the kinetic regime is independent of runoff in the generalized model. Moreover, the impact of a linear T -dependent runoff on thermodynamic weathering at a constant CO 2 partial pressure ( Figure 9a) is significant but not strong enough to mitigate the exponential effect of T on thermodynamic concentrations. Thus, the T -dependence of runoff, unless stronger than a linear dependence, does not impact the understanding of the T -sensitivity of weathering.

Global Silicate Weathering Rates
During the Archean (2.5−4 Ga), the incident solar radiation was about 70−80% of its present-day value, not high enough to maintain a temperate climate with present-day atmospheric CO 2 levels (Sagan & Mullen 1972;Charnay et al. 2020). Although there are no direct measurements of historical weathering rates, the Archean geologic record suggests a steady decrease in P CO2 from the order of 0.1 bar to modern values while maintaining surface temperatures between 280 K and 315 K (Krissansen-Totton et al. 2018, and references therein). The lower insolation during the Archean was compensated by the greenhouse effect of CO 2 . As the insolation increased, climates should have become warmer than the observations suggest. Without a negative feedback of silicate weathering that allowed a decrease in CO 2 levels as insolation increased, modern climates would not have been temperate (Walker et al. 1981;Berner et al. 1983;Kasting et al. 1993). The extent of silicate weathering during the history of Earth and the contribution of seafloor weathering is debated (e.g., Sleep & Zahnle 2001;Berner & Kothavala 2001;Foley 2015;Krissansen-Totton et al. 2018;Coogan & Gillis 2018). By applying the generalized weathering model to Earth, we lay the foundation for understanding climate regulation on temperate planets.
Since present-day continents on Earth are mostly granitic and the seafloor is basaltic, we choose a granite lithology for continents and a basalt lithology for the seafloor. Continental and seafloor weathering rates on Earth are calculated up to 4 Ga (Figure 10c,d). This is the first application of a transport-controlled weathering model to seafloor weathering on Earth. Since there are no direct measurements of historical weathering rates on Earth, we compare the generalized model developed in this study with a model from Krissansen-Totton et al. (2018, Figure 3, hereafter, KT18). Instead of using a climate model, time-dependent median models of P CO2 and T from the same KT18 model are used as inputs to the generalized model. The same inputs are used to obtain the continental weathering rate from Walker et al. (1981) and the seafloor weathering from Brady & Gíslason (1997), where the present-day weathering rates of both models are normalized to those of KT18 models.
Ideally, continental and seafloor weathering rates depend on the continental surface area (f A) and the seafloor surface area ((1 − f ) A), respectively (Equation 4). However, not all continental and seafloor surface area undergoes weathering on Earth. Fekete et al. (2002) report the value of continental weatherable area to be equal to 93 Mm 2 , about 60% of modern continental area (= 0.6 f A). About 147 Mm 2 of the seafloor area (= 0.41 (1 − f ) A) is expected to contribute to seafloor weathering (given by the exposed area of low-temperature hydrothermal systems, Johnson & Pruis 2003). Therefore, in Figure 10, the continental and seafloor weathering rates are given by W cont = w cont × 93 Mm 2 and W seaf = w seaf × 147 Mm 2 , respectively. A more precise definition of continental weatherable area is the area susceptible to precipitation and runoff, and that of seafloor weatherable area is the area covered by low-temperature hydrothermal systems that are younger than approximately 60 Myr (Stein & Stein 1994;Johnson & Pruis 2003;Coogan & Gillis 2018). Nonetheless, the definition given in Equation (4) is the best proxy for calculating global weathering rates on exoplanets, essentially giving upper estimates.
The present-day continental weathering rate for young soils (t s = 0) in Figure 10(c) is at the transition between thermodynamic and kinetic regimes (see also Figure 8a). Similarly, the present-day seafloor weathering rate for young pore-space is in the thermodynamic regime (Figure 10d). As the CO 2 levels increase, weathering is driven by kinetics which is strongly sensitive to T (Figure 10a). The sensitivity of our granite and basalt models to P CO2 is weak because k eff of these rocks obtained from Palandri & Kharaka (2004) is almost constant as a function of pH in basic solutions resulting from our cal- (c) Generalized young soils (ts = 0) and old soils (ts = 100 kyr) granite models take the same values for ψ = ψ0, q = 0.3 m yr −1 (modern mean runoff) and continental weatherable area of 93 Mm 2 (Fekete et al. 2002). (d) Generalized young pore-space (ts = 0) and old pore-space (ts = 50 Myr) basalt models take different values for ψ (ψ = 100ψ0 and ψ = 15ψ0, respectively) but the same value for q = 0.05 m yr −1 (modern mean hydrothermal fluid flow rate) and seafloor weatherable area of 147 Mm 2 (Johnson & Pruis 2003). Gray shaded regions are 95% confidence intervals of KT18 models. Pink and light green vertical lines are uncertainties in the estimates of present-day silicate weathering rates given by geological measurements and models. Colored disks denote the transition between thermodynamic and kinetic regimes, whereas inverted triangles represent the transition between kinetic/supply and CO2 dissolution regimes.
culations. This P CO2 dependence is stronger or weaker depending on the choice of minerals to model the rock.
The Walker et al. (1981) continental weathering rate and the Brady & Gíslason (1997) seafloor weathering rate exhibit a monotonic rise because of a stronger dependence on P CO2 (β = 0.3 and β = 0.23, respectively) than T . The continental weathering rate for young soils and the seafloor weathering rate for young pore-space are respectively about an order of magnitude and up to four orders of magnitude higher than that of KT18 models (Figure 10c,d). This is because not entire exposed planetary surface area contains fresh rocks for weathering. In fact, on Earth's continents, orogeny (mountain building) exposes fresh rocks to the surface that are highly susceptible to weathering (low t s ), whereas the contribution of cratons (ancient continental crust, high t s ) to weathering is smaller than mountains (Maher & Chamberlain 2014). Similarly, ridge volcanism exposes fresh basalt at midocean ridges but the majority of seafloor area is older than a million years. To match the present-day continental rate of our model with that of KT18, the soil age needs to be set to 100 kyr which is the modern mean soil age. To match the present-day seafloor weathering rate of our model with that of KT18, the pore-space age is set to 50 Myr (mean seafloor age) and the dimensionless pore-space parameter is decreased from 100ψ 0 to 15ψ 0 . The pore-space parameter has never been discussed in the context of seafloor weathering. On continents, this parameter varies with topography and its mean value for present-day continents is not well-constrained as previous studies adopt values from 0.1−10 × that of the reference value used in this study (Maher 2010(Maher , 2011Maher & Chamberlain 2014;Winnick & Maher 2018).
In addition to pore-space parameters, the present-day regional variation in continental runoff (0.01−3 m yr −1 , Gaillardet et al. 1999;Fekete et al. 2002) and seafloor fluid flow rate (0.001−0.7 m yr −1 , Stein & Stein 1994;Johnson & Pruis 2003;Hasterok 2013) is more than two orders of magnitude. At fluid flow rates lower than the mean values (arid climates or low hydrothermal heat flux, see Figure 7), the models might escape kinetic/supply regimes and enter the thermodynamic regime. These models rely on global mean values t s , ψ and q to be good proxies for calculating global weathering rates, that is not necessarily a good approximation. For example, mountains contribute an order of magnitude more weathering flux than cratons (Maher & Chamberlain 2014). Additionally, basaltic regions on continents (omitted in our calculations) may contribute a weathering flux (per unit area) higher than that of granitic regions by a factor of 5 (Ibarra et al. 2016). An appropriate way to calculate global weathering rates on Earth is the summation of local weathering rates. This is an excellent strategy for Earth but not for exoplanets whose global pore-space properties and fluid flow rates are not measurable let alone regional properties. Hence, it is best to study exoplanets by applying the generalized weathering model on the basis of globally-averaged values for these parameters.
Holding pore-space parameters constant throughout Earth's dynamic history is another critical approximation. The uplift of the Himalayan-Tibetan Plateau in the past 40 million years has decreased the mean soil age which has likely increased the global weathering rates (Kump et al. 2000). The lithology of continents has also evolved over time and therefore application of one type of lithology may result in inaccurate weathering rates for certain periods of Earth's history. The continental area fraction (and consequently the weatherable area) was smaller in the Archean than today. Even if the strength of continental weathering flux was similar at t = 4 Ga compared to today but the continental area was significantly smaller, then the continental weathering rate must have been significantly smaller. The continental and seafloor DIC fluxes for the old pore-space models enter the CO 2 dissolution regimes at ∼1.5 Ga and ∼0.7 Ga, respectively where the aqueous CO 2 concentrations dominate the total DIC (Figure 10c,d). However, the true contribution of [CO 2 (aq)] to silicate weathering depends on the carbonate compensation depth in oceans (Pytkowicz 1970) as well as reverse weathering (Mackenzie & Garrels 1966;Isson & Planavsky 2018; Krissansen-Totton & Catling 2020) that requires knowledge of ocean alkalinity and ocean pH, which suggests that an ocean chemistry model should form the basis of future research.

Weathering Regimes and the Habitable Zone
By measuring the amount of CO 2 in the atmospheres of Earth-size planets, it may be possible to statistically differentiate between runaway greenhouse (Venuslike) and temperate (Earth-like) planets on both sides of the inner-edge of the habitable zone (Bean et al. 2017;Checlair et al. 2019;Graham & Pierrehumbert 2020). The inner boundary of the habitable zone is defined as the distance from the host star, where the surface temperature of a planet increases to the critical point of water. Beyond this point, the presence of liquid water is, therefore, physically impossible. Reaching such a temperature can either be the result of a water vapor feedback mechanism, known as the runaway greenhouse, or by an additional greenhouse gas that pushes the surface temperature beyond the critical point of water on its own (the Simpson-Nakajima and Komabayasi-Ingersoll limits, Simpson 1929;Komabayasi 1967;In-gersoll 1969;Nakajima et al. 1992). Nakajima et al. (1992); Abe (1993) show that an additional greenhouse gas (e.g., CO 2 ) alone does not impact the critical insolation where the water-induced runaway greenhouse effect occurs. Only if that additional greenhouse gas alone can push the surface temperature of the planet beyond the critical point of water, would it impact on the location of the inner boundary of the habitable zone. If most liquid water ends up in the gaseous phase before the runaway greenhouse, there will be no such impact.
Assuming that silicate weathering feedback influences the extent of the habitable zone, the negative feedback of silicate weathering would move the inner-edge of the habitable zone towards the star, thereby extending its width. Normally, kinetics of weathering reactions is attributed to this negative feedback ( Figure 11). However, this feedback may not necessarily be negative. For instance, if weathering is limited by the supply of fresh rocks, the negative feedback is lost (West 2012;Foley 2015). In the supply regime, instead of climate or transport parameters, weathering depends on pore-space parameters that are held constant. However, if there exists a positive correlation between pore-space parameters and T , the negative feedback should be active even in the supply regime. In the thermodynamic regime, the weathering flux decreases with T and increases with x CO2(g) . If the T sensitivity dominates over the x CO2(g) sensitivity, the resulting thermodynamic weathering flux decreases with T . In this case, silicate weathering can provide positive feedback to the carbon cycle. Such a positive feedback may shift the inner-edge of the habitable zone away from the star (Figure 11). It is essential to elucidate the implications of a possible positive feedback on the onset of runaway greenhouse. Moreover, since the surface lithology strongly influences the weathering regimes, lithology may also impact the extent of the habitable zone. Fur-thermore, the role of weathering in predicting climate is intertwined with the assumption of weathering in climate models that provide a relation between CO 2 partial pressure, surface temperature, incident stellar flux, planetary albedo and the planet's ability to emit thermal radiation (Kasting et al. 1993;Kopparapu et al. 2013). The climate model calculations need to be revisited with a generalized weathering model based on lithology such as implemented in this study.

SUMMARY AND CONCLUSIONS
Silicate weathering is a key process in the carbon cycle that transfers CO 2 from the atmosphere to the surface of a planet. The intensity of silicate weathering has previously been attributed to kinetics of fluid-rock reactions (e.g., Walker et al. 1981;Berner et al. 1983;Kump et al. 2000). Maher (2011);Maher & Chamberlain (2014) show that if the reaction rate exceeds the fluid flow rate, thermodynamics of fluid-rock reactions at chemical equilibrium drives weathering instead of kinetics. This transport-controlled approach models both thermodynamic and kinetic regimes of weathering with a single formulation. Moreover, if there is a limited supply of fresh rocks, the weathering is supply-limited. The applications of this approach to continental weathering on Earth (Winnick & Maher 2018) and temperate exoplanets (Graham & Pierrehumbert 2020) consider weathering reactions of individual minerals.
In this study, we extend this approach to the weathering of any rock type (lithology) and apply it to seafloor weathering in addition to continental weathering. We demonstrate that it is essential to simultaneously consider weathering reactions of all minerals present in a rock as well as the reactions in the water-bicarbonate system instead of weathering reactions of individual minerals. Moreover, this model allows the calculation of absolute weathering rates instead of weathering rates normalized to present-day values as most previous weathering studies. In addition to climate properties (T and x CO2(g) ) and runoff or fluid flow rate (q), this model is mainly sensitive to age of soils (t s ) and a dimensionless scaling parameter (ψ) based on pore-space and rock properties. The equilibrium constants and kinetic rate coefficients are a function of T and x CO2(g) . Depending on these five parameters, the weathering for a given lithology is in the thermodynamic, kinetic, supply or CO 2 dissolution regimes. Close to the regime transition points, the contribution of both regimes to weathering is comparable.
Weathering reactions at chemical equilibrium give the maximum concentrations of weathering products. This approach is used to calculate an upper limit to weath-ering flux for a given lithology. The larger the fraction of divalent cations in rocks, the higher the sensitivity of maximum weathering flux to CO 2 volume mixing ratio. This thermodynamic x CO2(g) sensitivity (power-law exponent β) is 0.71, 0.65 and 0.54 for peridotite, basalt and granite, respectively. These values are subject to change depending on the choice of minerals to define a rock type. The thermodynamic x CO2(g) sensitivity of these rocks is stronger than the kinetic x CO2(g) sensitivity implemented in previous studies (0.22−0.55, Walker et al. 1981;Berner 1991;Driscoll & Bercovici 2013). However, the combined effect of x CO2(g) and T results in a stronger kinetic x CO2(g) sensitivity. The thermodynamic T sensitivity of weathering of rocks is negative, implying that weathering decreases with T . This is in contrast with the prevalent understanding that weathering intensifies with an increase in T which is attributed to an increase in kinetic rate coefficients of mineral dissolution reactions with T (Lagache 1965;Palandri & Kharaka 2004;Brantley et al. 2008). An important implication of this finding is that when T increases without a strong variation in x CO2(g) , silicate weathering has the potential to instigate a positive feedback to the carbon cycle.
The transport-controlled model demonstrates that the weathering flux cannot necessarily be approximated by the kinetic weathering expression (Equation 2). Planets with mean arid climates (low runoff) and elevated topography (young soils) are likely in the thermodynamic regime of weathering. Moreover, limited supply of fresh rocks mitigates the role of kinetics. Application of this model to Earth suggests that global mean continental granite and seafloor basalt weathering rates are likely limited by the supply of fresh rocks, yet regional weathering rates can be influenced by both kinetics and thermodynamics. Studies modeling the extent of the habitable zone are based on the presumption that silicate weathering provides a negative feedback to the carbon cycle owing to kinetics of weathering reactions. This statement does not hold in the thermodynamic regime if T exerts a stronger control on weathering than x CO2(g) (positive feedback) or in the supply regime (no feedback). The focus of future studies should be on applying a generalized weathering model encompassing multiple weathering regimes to model climate as well as the boundaries of the habitable zone. The equilibrium constant K of a reaction is given by the difference of the Gibbs energy of formation of products and reactants as follows where ∆ f G P,T,i is the Gibbs energy of formation of i th species at pressure P and temperature T , ν i is the stoichiometric coefficient and R is the universal gas constant. The Gibbs energy of formation of each species is computed at any P and T in terms of the Gibbs energy of formation at reference pressure P 0 and reference temperature T 0 by where S P0,T0 is the entropy at the reference pressure and temperature, C P is the heat capacity at constant pressure as a function of temperature, and V is volume as a function of pressure. In this study, ∆ f G values are obtained from the CHNOSZ database (Dick 2019). The equilibrium constants of reactions given in Table B.1 are shown as a function of P and T in Figure A.1.
The kinetic rate coefficient k eff of mineral dissolution reactions are obtained from the compilation of Palandri & Kharaka (2004). This compilation (Table A. mental kinetics data, where E acid , E neut and E base are the activation energies at acidic, neutral and basic pH, A acid , A neut and A base are the preexponential factors at acidic, neutral and basic pH, and n acid and n base are the neutral and basic power-law exponents. Figure A.2 shows the variation of k eff with T and pH for a number of mineral dissolution reactions. For minerals not present in this compilation (ferrosilite, annite, grunerite), the kinetic rate coefficients are obtained from respective endmember minerals of the same mineral group (enstatite, phlogopite, anthophyllite).

B. MAXIMUM AND GENERALIZED CONCENTRATIONS OF ROCKS AND MINERALS
Sects. 2.3 and 2.4 introduce the methods to compute maximum (thermodynamic) and generalized solute concentrations for peridotite weathering. Table B.1 lists the mineral dissolution (reactions (a−n)) and waterbicarbonate reactions (reactions (o−q)), and the relation between equilibrium constants of these reactions and thermodynamic activities. Table B.2 gives the polynomial equations to calculate the activity of HCO − 3 at chemical equilibrium as a function of x CO2(g) for the weathering of all rocks and minerals considered. As an example, the maximum [HCO − 3 ] for peridotite weathering is obtained as a function of CO 2 volume mixing ratio, surface temperature and total pressure in Figure B.3. As described in Section 2.3, [HCO − 3 ] eq is sensitive to x CO2(g) and T . However, P has a negligible effect on [HCO − 3 ] eq because the equilibrium constants of reactions are largely unchanged up to 1000 bar (Appendix A). This figure demonstrates that the total pressure plays a negligible role in determining the solute concentrations of aqueous species. The effect of precipitation of amorphous silica (Winnick & Maher 2018) is not modeled since this effect changes the weathering flux by less than an order of magnitude (only at high x CO2(g) ) which is smaller than the 5−10 orders of magnitude spread in the weathering fluxes discussed in this study.
Once the maximum [HCO  6). This transport-controlled or diluted [HCO − 3 ] is then used as an input to solve for concentrations of other aqueous species such as CO 2− 3 , H + and OH − by assuming that the water-bicarbonate reactions obey chemical equilibrium. Figure B.2 demonstrates that generalized solute concentrations (Section 2.4) are strongly sensitive to lithology. For example, the transition between thermodynamic and kinetic weathering regimes of peridotite occurs at x CO2(g) = 10 ppmv for q = 0.3 m yr −1 . Whereas, this transition occurs at x CO2(g) > 10 5 ppmv for anorthite. Once these generalized concentrations are obtained, the generalized weathering flux is calculated using Equation (3)  Schematic describing the methodology of the weathering model CHILI. Square boxes represent software modules, ovals denote parameters (Tables 1 and 3) and rounded squares represent output quantities ( Table 2). The solute transport equation of Maher & Chamberlain (2014) is implemented to calculate transport-controlled solute concentrations. Thermodynamics and kinetics data are obtained from Dick (2019) and Palandri & Kharaka (2004), respectively (see Appendix A). Dissolved inorganic carbon components and pH as a function of x CO 2 (g) at T = 288 K (modern surface temperature), q = 0.3 m yr −1 (modern mean runoff) and ts = 0 (young soils). From left to right, colored disks mark the transition between thermodynamic and kinetic regimes, and the inverted triangle marks the transition between kinetic and CO2 dissolution regimes with respect to DIC. Reactions and relations between equilibrium constants and activities of reactants and products.
10 −2 10 0 10 2 10 4 10 6 x CO 2 (g) [ppmv] Henry's law states that the amount of gas dissolved in the liquid ([CO 2 (aq)] = a CO2(aq) × 1 mol dm −3 ) is proportional to its partial pressure above the liquid, P CO2 = x CO2(g) P , where x CO2(g) is the CO 2 volume mixing ratio and P is the total pressure. For the CO 2 dissolution reaction ((o) in Table B.1), the proportionality constant is the equilibrium constant K CO2 that itself depends on pressure and temperature, a CO2(aq),thermo = K CO2 x CO2(g) . (C7) We obtain the dimensionless K CO2 as a function of T and P from the CHNOSZ thermodynamic database (Dick 2019). Pierrehumbert (Equation 8.14, 2010) provide an approximate dimensional Arrhenius-type fitting function K H at any temperature T for Henry's law constant, with empirical factors, K 0 H = 1600 bar mol water mol CO2(aq) at a reference temperature T 0 = 298 K, and C H = 2400 K. The relation between a CO2(aq) and x CO2(g) using K H is given by a CO2(aq),arrhen = u K H x CO2(g) , where u = 55.5 bar mol water mol CO2(aq) (1 dm 3 of water contains 55.5 moles of water) is a conversion factor between the standard states of K H (1 bar and 1 mol CO2(aq) mol water ) and K CO2 (1 bar and 1 mol CO2(aq) dm 3 water ). In Figure C.1, these two models (Equations C7 and C9) are compared with the fit to experimental data on the solubility of CO 2 in pure water compiled by Diamond & Akinfiev (2003). The Arrhenius-type model is within 6% of that of the experimental data up to 330 K and deviates by up to 37% at higher temperatures. The thermodynamic model performs better than the Arrhenius-type model at all temperatures except for 288−313 K and is within 10% of the experimental data at 373 K. For this reason, we use the thermodynamic model to calculate the solubility of CO 2 in water instead of the Arrhenius-type model (Equation C8). Difference between model and experimental data for a CO 2 (aq) as a function of T . The experimental data is obtained from the compilation of Diamond & Akinfiev (2003). The models are given by the Arrhenius-type equation (Pierrehumbert 2010) and the thermodynamic database (Dick 2019).

D. SENSITIVITY OF THE DAMKÖHLER COEFFICIENT TO PARAMETERS
The Damköhler coefficient D w depends on seven parameters and two variables (Equation 7).
The two variables, equilibrium solute concentration C eq (= [HCO − 3 ] eq in this study) and effective kinetic rate coefficient k eff , are treated as free parameters in Figure D.1. All nine quantities are varied for a maximum possible range of their known values ( Figure D.1). D w is largely sensitive to four quantities, C eq , k eff , t s (age of soils) and L (flowpath length). The flowpath length is absorbed into the dimensionless pore-space parameter ψ which is a control parameter for the models in the main text (Equation 7). Figure D.1(b,c) highlights the interdependence of k eff and t s . At low k eff or low t s , D w is strongly sensitive to k eff and insensitive to t s , implying the presence of 'fast kinetic' regime. At high k eff or high t s , D w is independent of k eff and decreases strongly with t s , implying that weathering is in the 'slow kinetic' or supply-limited regime due to insufficient supply of fresh rocks for weathering. Figure D.1 also compares D w of our granite model to that of Kopparapu et al. (2014). These two models show similar trends between D w and respective parameters. The difference between the two models arise mainly from our assumption of endmember silicate minerals, instead of solid solutions, for the granite model that creates about an order of magnitude difference between C eq values.

E. CLIMATE MODELS
A climate model provides a relation between the surface temperature T , the CO 2 partial pressure P CO2 , top-of-atmosphere stellar flux S and planetary albedo α. Kadoya & Tajika (2019) provide a fitting function to the climate model of Kopparapu et al. (2013Kopparapu et al. ( , 2014 which is valid for T in the range 150−350 K with P CO2 in the range 10 −5 − 10 bar at saturated H 2 O and 1 bar N 2 . The fitting function is given by T = [1 ξ ξ 2 ξ 3 ξ 4 ξ 5 ξ 6 ], (E11) where the outgoing longwave radiation F OLR is a function of T and P CO2 , I 0 = −3.1 W m −2 and ξ = 0.01 (T − 250). For P CO2 < 1 bar, For P CO2 > 1 bar, This fit assumes that the total pressure is the sum of partial pressures of CO 2 , N 2 and H 2 O, P = P CO2 + P N2 + P H2O . Then the volume mixing ratio of CO 2 is given by x CO2(g) = PCO 2 PCO 2 +PN 2 +PH 2 O . 2469 J g −1 is the latent heat of vaporization for water at the standard temperature and R = 8.3145 J K −1 mol −1 is the universal gas constant. Equations (E10−E16) are solved by balancing the energy fluxes of the globally-averaged absorbed instellation S avg and F OLR , where and with the geometric factor 4 comes from ratio of the planet surface area to the area of its cross-section. For present-day albedo (α = 0.3) and present-day solar flux (S = 1360 W m −2 ), this fit results in T between 280 K and 350 K and P CO2 between 10 −5 bar and 0.5 bar (equivalently x CO2(g) between 10 −5 − 0.2). Another climate model used frequently in carbon cycle studies (e.g., Foley 2015) is the one from Walker et al. (1981). The relation between T and P CO2 is given by T = T * + 2(T e − T * e ) + 4.6 P CO2 P * where T e is the effective temperature given by T e = (S avg /σ SB ) 1/4 with σ SB = 5.67 × 10 −8 W m −2 K −4 as the Stefan-Boltmann constant and the present-day values of temperature, effective temperature and CO 2 partial pressure are assumed to be T * = 285 K, T * e = 254 K and P * CO2 = 330 × 10 −6 bar, respectively (Kasting et al. 1984). Figure E.1 shows the comparison between the models of Kadoya & Tajika (2019) and Walker et al. (1981). Both models show almost the same temperatures for α between 0.3 and 0.5 at P CO2 = 280×10 −6 bar. However, for α < 0.25, the Kadoya & Tajika (2019) model shows a steep temperature rise with decreasing α. As a function of P CO2 at α = 0.3, both models exhibit temperatures within 5% of each other.