Babak Ravaji

Víctor Alí-Lagoa

Marco Delbo

J W Wilkerson

J Mike Walker '

Unraveling the mechanics of thermal stress weathering: rate-effects, size-effects, and scaling laws

Keywords: Geological processes, Asteroids, Moon, Thermal stress, Space weathering, Multi-scale modeling

Thermal stress weathering is now recognized to be an active and significant geomorphological process on airless bodies. This study aims to understand the key factors governing thermal stresses in rocks on airless bodies through extensive numerical calculations and analytic analyses. Some of the key factors governing thermal stresses, are found to be the diurnal surface temperature variation, the second-order spatial gradient of the temperature field, the thermal skin depth, and the rock size of interest. Microscopic (grain-scale) thermal stresses are driven primarily by the amplitude of the magnitude of the maximum diurnal temperature variation at said depth. Macroscopic (rock-scale) thermal stresses are more complex, and their nature fundamentally depends on two length scales: the thermal skin depth and the rock size of interests. For rock sizes larger than the thermal skin depth, macroscopic thermal stresses are driven primarily by second (and higher) order spatial gradients of temperature. For rock sizes smaller than the thermal skin depth, macroscopic thermal stresses are primarily driven by the ratio of rock size to thermal skin depth with macroscopic thermal stresses being greatest when this ratio is 1/2. Additionally, scaling relations for diurnal surface temperature variation, time-rateof-change of surface temperature, as well as peak microscopic (grain-scale) and macroscopic (rock-scale) thermal stresses are derived to provide a more accessible modeling tool. These scaling relations are remarkably accurate when compared to both the numerical calculations as well as three-dimensional finite element calculations. The model formulation, results, and scaling relations provided here allow the estimation of diurnal temperatures and thermal stresses on rocks of various size and materials on airless bodies at any orbital distance with a broad spectrum of spin rates. Lastly, we postulate and confirm that there is a critical spin rate where macroscopic thermal stresses will be greatest.

Introduction

In most terrestrial environments, geomorphology is governed by a combination of physical, chemical, and biological processes with aeolian, fluvial, frost, and chemical processes being the primary drivers [START_REF] Wellman | Salt Weathering, a Neglected Geological Erosive Agent in Coastal and Arid Environments[END_REF][START_REF] Cooke | Salt Weathering in Deserts[END_REF][START_REF] Hallet | Why Do Freezing Rocks Break? Science[END_REF][START_REF] Veverka | Planetary geology in the 1980s[END_REF]. In arid terrestrial environments many of these processes are inactive [START_REF] Viles | Stress histories control rockbreakdown trajectories in arid environments[END_REF][START_REF] Hall | The role of thermal stress fatigue in the breakdown of rock in cold regions[END_REF], and even fewer are active in airless, extraterrestrial environments. Surface morphology evolution on airless planetary bodies, such as the Moon, Mercury, and asteroids, has historically been primarily attributed to impact cratering, micrometeorite bombardment [START_REF] Housen | Asteroidal regoliths[END_REF][START_REF] Horz | Impact experiments related to the evolution of planetary regoliths[END_REF]; Collisions with solar and galactic cosmic rays, irradiation, ˚Corresponding author

Email address: wilkerson@tamu.edu (J.W. Wilkerson) implantation, and sputtering from solar wind particles also cause a space weathering of these surfaces [START_REF] Brunetto | Asteroid Surface Alteration by Space Weathering Processes[END_REF]. Despite the fact that these airless bodies experience large diurnal temperature variations (on the order of " 100 K), thermal stress weathering has long been presumed to be of little significance in the inner solar system. Within the last decade, thermal stress weathering has been revisited with greater rigour and is now suspected to play an important role in rock breakdown, regolith generation, crater degradation, and landscape evolution in Earth's deserts and cold regions [START_REF] Hall | The role of thermal stress fatigue in the breakdown of rock in cold regions[END_REF][START_REF] Lamp | Thermal stress weathering and the spalling of Antarctic rocks[END_REF], Mars [START_REF] Viles | Simulating weathering of basalt on Mars and Earth by thermal cycling[END_REF][START_REF] Eppes | Cracks in Martian boulders exhibit preferred orientations that point to solar-induced thermal stress[END_REF], Mercury (Molaro and Byrne, 2012), Moon [START_REF] Molaro | Grain-scale thermoelastic stresses and spatiotemporal temperature gradients on airless bodies, implications for rock breakdown[END_REF][START_REF] Mazrouei | Does Thermal Fatigue Play a Role in Lunar Regolith Formation[END_REF][START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF], near-Earth asteroids [START_REF] Graves | Resurfacing asteroids from thermally induced surface degradation[END_REF][START_REF] Delbo | Thermal fatigue as the origin of regolith on small asteroids[END_REF][START_REF] Andrew J Dombard | Boulders and ponds on the Asteroid 433 Eros[END_REF][START_REF] Jewitt | The active asteroids[END_REF], and perhaps comets [START_REF] Shestakova | The thermal destruction of solids near the Sun[END_REF][START_REF] L V Tambovtseva | Non-LTE models of the accretion disks of UX Ori stars[END_REF][START_REF] Pajola | The pristine interior of comet 67P revealed by the combined Aswan outburst and cliff collapse[END_REF][START_REF] Ramy El-Maarry | Regional surface morphology of comet 67P/Churyumov-Gerasimenko from Rosetta/OSIRIS images[END_REF][START_REF] V Alí-Lagoa | Rapid temperature changes and the early activity on comet 67P/ChuryumovGerasimenko[END_REF]. Field observations of thermal stress weathering on Earth have been well-documented (Ollier, 1969;[START_REF] Rice | Insolation warmed over[END_REF][START_REF] Siegesmund | Physical weathering of marbles caused by anisotropic thermal expansion[END_REF][START_REF] Koch | The combined effect of moisture and temperature on the anomalous expansion behaviour of marble[END_REF][START_REF] Weiss | Insolation weathering and hygric dilatation: two competitive factors in stone degradation[END_REF][START_REF] Heather | Microclimate and weathering in the central Namib Desert, Namibia[END_REF][START_REF] Sumner | Rock surface temperatures in southern Namibia and implications for thermally-driven physical weathering[END_REF]Eppes et al., 2010;[START_REF] Aldred | The influence of solar-induced thermal stresses on the mechanical weathering of rocks in humid mid-latitudes[END_REF][START_REF] Brian | Rockfall triggering by cyclic thermal stressing of exfoliation fractures[END_REF][START_REF] Lamp | Thermal stress weathering and the spalling of Antarctic rocks[END_REF]. In particular, [START_REF] Mcfadden | Physical weathering in arid landscapes due to diurnal variation in the direction of solar heating[END_REF], Eppes et al. (2010) and (Eppes et al., 2016) observed that most crack features in boulders lying in the deserts of Arizona, east United States, Australia, and Mongolia exhibited an N-S orientation, strong circumstantial evidence that these crack features were driven by thermal stresses induced by heating from the sun moving along an E-W path relative to the boulder. Similar crack feature orientations have be observed by the rovers on Mars [START_REF] Eppes | Cracks in Martian boulders exhibit preferred orientations that point to solar-induced thermal stress[END_REF]. If thermal stress weathering is operative in bodies with atmospheres, then it is very likely more effective on airless bodies where the lack of an atmosphere enables larger diurnal temperature variations and faster rate of temperature change [START_REF] Jewitt | Activity in Geminid Parent (3200) Phaethon[END_REF]. However, even in arid terrestrial environments there remains some degree of doubt regarding the assumption that thermal stress weathering is dominate over the freeze-thaw mechanism [START_REF] Carpinteri | Are the Paris Law parameters dependent on each other?[END_REF][START_REF] Hall | The role of thermal stress fatigue in the breakdown of rock in cold regions[END_REF]. This assumption has been addressed, to some degree, through laboratory experiments. The earliest primitive laboratory experiments were carried out by [START_REF] David | The factor of fatigue in rock exfoliation[END_REF] and [START_REF] Blackwelder | The insolation hyphothesis of rock weathering[END_REF], who concluded that thermal stress weathering is likely insignificant, unless rapid cooling is achieved via water quench. However, more modern laboratory experiments of [START_REF] Luque | Direct observation of microcrack development in marble caused by thermal weathering[END_REF][START_REF] Levi | Thermal fatigue: a possible source of structural modifications in meteorites[END_REF][START_REF] Delbo | Thermal fatigue as the origin of regolith on small asteroids[END_REF][START_REF] Brian | Rockfall triggering by cyclic thermal stressing of exfoliation fractures[END_REF][START_REF] Hazeli | The origins of Asteroidal rock disaggregation: Interplay of thermal fatigue and microstructure[END_REF] have demonstrated that thermal stress weathering is operative even in anhydrous environments.

Most of our classical understanding of rock breakdown and subsequent regolith evolution on airless bodies in the solar system has been gleaned from observations on the Moon coupled with studies of returned samples of lunar regolith [START_REF] Gault | Effects of microcratering on the lunar surface[END_REF][START_REF] Gault | Experimental studies of oblique impact[END_REF][START_REF] Basilevsky | Survival times of metersized boulders on the surface of the Moon[END_REF][START_REF] Krishna | Impact spallation processes on the Moon: A case study from the size and shape analysis of ejecta boulders and secondary craters of Censorinus crater[END_REF]. However, recent observations on near-Earth asteroids (NEAs), e.g. Itokawa [START_REF] Yano | Touchdown of the Hayabusa spacecraft at the Muses sea on Itokawa[END_REF] and 433 Eros [START_REF] Bell | Near-IR Reflectance Spectroscopy of 433 Eros from the NIS Instrument on the NEAR Mission[END_REF][START_REF] Clark | Space weathering on Eros: Constraints from albedo and spectral measurements of Psyche crater[END_REF][START_REF] Murchie | Color Variations on Eros from NEAR Multispectral Imaging[END_REF], and main-belt asteroids (MBAs), e.g. 243 Ida [START_REF] Veverka | Head III. Ida and Dactyl: Spectral Reflectance and Color Variations[END_REF][START_REF] Clark | S-type asteroids, ordinary chondrites, and space weathering: The evidence from Galileo's fly-bys of Gaspra and Ida[END_REF], suggest that regolith generation and evolution on asteroids is quite different from that of the Moon [START_REF] Michael | Space weathering and the interpretation of asteroid reflectance spectra[END_REF]. Therefore, in most of the cases expanding our Moon understandings and observations to other solar system objects, e.g. NEAs, is not necessarily appropriate.

In the case of asteroids, tests and observations directly on these bodies are extremely valuable, but extremely rare. Numerical and analytical studies in conjunction with remote experiments in this area would increase our understanding of their surface evolution. In order to assess the effects of thermal cycles on the surface of asteroids, experiments have been conducted so far on asteroid analogs. [START_REF] Levi | Thermal fatigue: a possible source of structural modifications in meteorites[END_REF] confirmed the thermal fatigue hypothesis in specimens of two H-chondrites. Recently Molaro and Byrne (2012); [START_REF] Molaro | Grain-scale thermoelastic stresses and spatiotemporal temperature gradients on airless bodies, implications for rock breakdown[END_REF][START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF] performed series of finite element studies to show the significance of thermal stress on airless bodies. [START_REF] Delbo | Thermal fatigue as the origin of regolith on small asteroids[END_REF]; [START_REF] Hazeli | The origins of Asteroidal rock disaggregation: Interplay of thermal fatigue and microstructure[END_REF] experimentally tested and conducted numerical studies to respectively determine thermal cycling degradation rates of meteorites and regolith on NEAs. These works suggest that the lifetime of surface rocks against thermal fatigue is, in general, shorter than the lifetime due to breakdown by micrometeorite bombardment.

However, recently [START_REF] Basilevsky | Survival times of meter-sized rock boulders on the surface of airless bodies[END_REF] raised doubts regarding the effectiveness of thermal fatigue. They note that lunar observations show large boulders with no clear signs of thermal degradation even after 20 ´400 M a. [START_REF] Basilevsky | Survival times of meter-sized rock boulders on the surface of airless bodies[END_REF] conclude that this is evidence of irregularities in the model of [START_REF] Delbo | Thermal fatigue as the origin of regolith on small asteroids[END_REF], which calculates thermal degradation rates on fast rotating (6 hour periods) NEAs. In their claim, [START_REF] Basilevsky | Survival times of meter-sized rock boulders on the surface of airless bodies[END_REF] approximate thermal stresses as being independent of the rotation period, and thus assume that the degradation rates per thermal cycle on fast rotating (6 hour periods) bodies is approximately the same as on very slow rotating bodies with "700 hour periods. This was shown to be poor approximation by [START_REF] Ravaji | The effect of rotation period on thermal stress weathering[END_REF], who demonstrated that macroscopic thermal stresses on a body with a "700 hour period may be an order of magnitude smaller than on a body with a 6 hour period. Thus, the lunar observations noted by [START_REF] Basilevsky | Survival times of meter-sized rock boulders on the surface of airless bodies[END_REF] are not necessarily evidence of irregularities in the model of [START_REF] Delbo | Thermal fatigue as the origin of regolith on small asteroids[END_REF], since the macroscopic thermal stresses are so different in the two cases. Here we derive a more accurate scaling approach that accounts for the strong dependence of spin rate on thermal stresses.

The efficacy of thermal stress weathering has traditionally been believed to be driven largely by two main factors: (i) diurnal surface temperature variation ∆T s [START_REF] Patrick | Effect of Cracks on Elastic Properties of Low Porosity Rocks[END_REF][START_REF] Yong | Thermally induced acoustic emission in Westerly granite[END_REF][START_REF] Molaro | Grain-scale thermoelastic stresses and spatiotemporal temperature gradients on airless bodies, implications for rock breakdown[END_REF] and (ii) rate of surface temperature change dT s {dt [START_REF] Patrick | Effect of Cracks on Elastic Properties of Low Porosity Rocks[END_REF][START_REF] Richter | Thermal expansion behavior of igneous rocks[END_REF][START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF]. Moreover, these two driving factors are themselves largely governed by many other conditions, e.g. orbital distance, rotation rate, orientation relative to the Sun, latitude, presence of shadows, atmosphere, and thermophysical material properties. For example, [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF] noted that thermal stress weathering is anticipated to be the most effective on slow rotating bodies, e.g. the Moon, as well as bodies with small orbital distances, e.g. Mercury and the NEA (3200) Phaethon which have perihelion distances of 0.3 and 0.14 AU, respectively. Here, we propose quantitative classifications of slow/fast rotating bodies from a thermal stress standpoint and derive simple scaling relations for thermal stresses that demonstrate strong agreement with the full three-dimensional finite element calculations of [START_REF] Molaro | Grain-scale thermoelastic stresses and spatiotemporal temperature gradients on airless bodies, implications for rock breakdown[END_REF][START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF].

The current paper aims to provide a more general understanding of thermal stresses in rocks on the surface of celestial bodies with a vast range of spins and orbital distances as well as thermophysical and thermoelastic properties. The continuum mechanics perspective on the thermal stresses is discussed in section 2. The thermal stress weathering model utilized in this paper is briefly outlined in section 3. Key factors governing diurnal temperature variations, time-rateof-change, spatial gradients, and temperature fields are discussed in section 4. In section 5 we explore implications of the key aspects of the spatiotemporal temperatures fields on thermal stresses. This section corroborates, builds upon, and attempts to generalize recent work of [START_REF] Molaro | Grain-scale thermoelastic stresses and spatiotemporal temperature gradients on airless bodies, implications for rock breakdown[END_REF][START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF], who first reported several of the key factors discussed here. In particular, our discussion and analytic analysis of macroscopic thermal stresses section 5 provides additional interpretations and clarity. Simple, yet accurate and powerful scaling relations provided in this sections are indeed fruits of such a deep fundamental approach. Finally, summary and conclusions are provided in section 6.

The continuum mechanics perspective on thermal stresses and terminology

In much of the literature on thermal stress weathering, the two main driving factors, i.e. ∆T s and dT s {dt, have been regarded as different mechanisms, typically termed thermal fracture/fatigue and thermal shock, respectively. However, these two mechanisms are, in fact, quite similar from a continuum mechanics standpoint. Both mechanisms involve the development of internal stresses that are required to ensure that the compatibility equation, c.f. (Sadd, 2014), of continuum mechanics is satisfied, i.e.

∇ ˆp∇ ˆεq " 0.

(1)

These compatibility equations pertain to the total strain tensor ε, which may be additively decomposed into elastic and thermal parts, respectively, i.e. ε " ε e `εt . The thermal strain tensor mathematically describes the thermal expansion field of the body, which for isotropic materials is simply ε t " α∆T I, where α is the local linear thermal expansion coefficient and I denotes the second-order identity tensor. The elastic strain tensor ε e describes the part of the deformation that is due to internal stresses (e.g. stretching resulting from an applied force). The total strain tensor relates to the actual observed deformation of the body, which is a combination of elastic and thermal contributions. Consider now a stress-free (isotropic) body, i.e. ε e " 0, that is free to expand in all directions (unconstrained). As such, the total strain tensor is exactly equal to the thermal strain tensor, i.e. ε " ε t " α∆T I. It is easily seen that both homogeneous and linear spatial variations of the thermal expansion field α∆T satisfy Eq. ( 1). However, non-linear spatial variations of α∆T do not generally satisfy this compatibility equation, implying that a non-zero elastic strain field, i.e. ε e ‰ 0, is required to satisfy Eq. ( 1) and hence a non-zero stress field. This stress field required to maintain compatibility is the so-called thermal stress.

The term thermal fracture/fatigue is intended to denote cases for which the thermal stresses are primarily produced by non-linear spatial variations in the local thermal expansion coefficient α. Typically, these spatial variations are on the grain-scale of the material where neighboring grains may have vastly different thermal expansion coefficients (and at the grain-scale the thermal expansion coefficients may be anisotropic). For a given material, the magnitude of these thermal stresses scales with the diurnal temperature variation, ∆T s , c.f. [START_REF] Lawson | Geodynamics[END_REF].

On the other hand, the term thermal shock is meant to convey that the thermal stresses are produced due to rapid changes in the surface temperature, i.e. |dT s {dt| " 0. Experimental observations have led much of the community to view |dT s {dt| ą 2 K/min as a critical threshold, above which rocks are susceptible to breakdown via thermal shock [START_REF] Patrick | Effect of Cracks on Elastic Properties of Low Porosity Rocks[END_REF][START_REF] Richter | Thermal expansion behavior of igneous rocks[END_REF]. However, the threshold of 2 K/min is criticized by more recent researchers, e.g. in experimental and numerical study conducted by Eppes et al. (2016), and emphasis is placed on the importance of spatial temperature gradient instead of surface temperature rate [START_REF] Boelhouwers | Critical assessment of the 2cmin1 threshold for thermal stress weathering[END_REF][START_REF] Molaro | Stress, on the rocks: Thermally induced stresses in rocks and microstructures on airless bodies, implications for breakdown[END_REF]. Indeed, it should be recognized that although dT s {dt is empirically related to the magnitude of thermal stress and the efficacy of thermal shock, it is not the fundamental governing factor. According to the heat conduction equation, i.e. ρc p dT dt " k∇ 2 T, with k, ρ and c p denoting thermal conductivity, density and heat capacity, respectively, a rapid surface temperature change will inevitably induce a spatially non-linear temperature field, and the magnitude of dT s {dt is (nearly) proportional to the magnitude of the second-order (and higher) spatial gradients of temperature. According to compatibility equation, i.e. Eq. ( 1), it is in fact, these second-order ∇ 2 T (and higher) spatial gradients of temperature that govern the thermal stresses associated with socalled thermal shock. This seems to have been recognized by [START_REF] Lu | The thermal shock resistance of solids[END_REF], but seems to have been misunderstood in some of the more recent literature. This distinction is particularly important when considering thermal shock of different materials. For example, the thermal diffusivity, k{ρc p , of metals is roughly " 100 times larger than that of rocks, so the same surface temperature heating rate of 2 K/min would generate thermal gradients, ∇ 2 T s , that are " 100 times smaller in metals as compared to rocks. This is fundamentally why a 2 K/min heating rate may crack rocks, but may be inconsequential for metals. That said, dT s {dt is obviously easier to measure in the laboratory, field, or through remote observations. Essentially all fracture and fatigue processes are activated by internal stresses within the body exceeding some material-dependent threshold measure. In order to satisfy compatibility, i.e. Eq. ( 1), thermal stress fields are necessarily induced by second-order (and higher) spatial gradients of α∆T . Depending on how those spatial gradients come about, the community has adopted different terminology, namely thermal fracture/fatigue (for ∇ 2 α ‰ 0 and ∇ 2 T " 0) and thermal shock (for |dT {dt| " pk{ρc p q ˇˇ∇ 2 T ˇˇÁ 2 K/min); however, this may have led to unnecessary confusion, as noted by [START_REF] Molaro | Grain-scale thermoelastic stresses and spatiotemporal temperature gradients on airless bodies, implications for rock breakdown[END_REF]. [START_REF] Delbo | Thermal fatigue as the origin of regolith on small asteroids[END_REF] have attempted to resolve this confusion by casting their framework in terms of microscopic (grain-scale) thermal stresses and macroscopic (boulder-scale) thermal stresses, both of which may operate simultaneously on airless planetary bodies in the inner solar system. Indeed, the macroscopic (rock-scale) stress is the driving force for N-S oriented cracks observed on the Earth surface [START_REF] Mcfadden | Physical weathering in arid landscapes due to diurnal variation in the direction of solar heating[END_REF]Eppes et al., 2010Eppes et al., , 2016)), while the microscopic one is the root of granular disintegration in Earth rocks [START_REF] Hall | Rock thermal data at the grain scale: applicability to granular disintegration in cold environments[END_REF][START_REF] Gómez-Heras | Surface temperature differences between minerals in crystalline rocks: Implications for granular disaggregation of granites through thermal fatigue[END_REF][START_REF] Cary | Granular disintegration of marble in nature: A thermal-mechanical origin for a grus and corestone landscape[END_REF]. The fact that both mechanisms may be operating simultaneously makes linear scaling of [START_REF] Basilevsky | Survival times of meter-sized rock boulders on the surface of airless bodies[END_REF] particularly problematic. A more rigorous scaling will be proposed in this paper.

A model for thermal stresses

For the present investigation, we make use of the thermal stress weathering model proposed by [START_REF] Delbo | Thermal fatigue as the origin of regolith on small asteroids[END_REF]. A brief overview of this model is provided here. The model considers an unconstrained spherical rock of diameter D on the surface of an asteroid, comet, moon, or planet with a spin rate of ω. The rock is considered part of a continuum regolith bed. As the body rotates in relation to the Sun, the rock experiences oscillating non-linear temperature fields, which generate thermal stresses throughout the rock. These thermal stresses, if of sufficient magnitude, can gradually (or catastrophically) drive crack growth.

In addition to these macroscopic geometric considerations, the model accounts for an idealized heterogeneous Fig. 1: Schematic of a spherical rock with an idealized microstructure composed of inclusions, e.g. chondrules, of radius r i positioned on a cubic lattice with nearest-neighbor spacing . At any macroscopic location z, there is a coresponding microscopic spatial coordinate y, whose origin lies at the center of the inclusion nearest to the macroscopic location z of interest.

microstructure composed of spherical inclusions, e.g. chondrules, of radius r i and spacing , as shown in Fig. 1. Here, we assume a microstructure similar to that of ordinary condrite materials, i.e. r i " 1.16 mm and " 1.4 mm according to [START_REF] Delbo | Thermal fatigue as the origin of regolith on small asteroids[END_REF].

Here, we limit the scope of the present paper to providing a more fundamental understanding of the thermal stresses that drive thermal stress weathering. The requisite equations of the [START_REF] Delbo | Thermal fatigue as the origin of regolith on small asteroids[END_REF] thermal stress weathering model are restated here for convenience.

Combined macroscopic and microscopic stresses:

σ " Σ `σT M , (2)
where σ denotes the equibiaxial in-plane pe c n direction) normal components of the total (combined) thermal stresses, which is additively decomposed into corresponding macroscopic pΣq and microscopic pσ T M q components. Σ denotes the macroscopic thermal stress associated with macroscopic non-linear thermal gradients, i.e. ∇ 2 T ‰ 0. σ T M denotes the corresponding normal component of the microscopic thermal stresses generated by heterogeneous thermal expansion coefficients in the microstructure, i.e. ∇ 2 α ‰ 0. Both the macroscopic and microscopic thermal stresses contribute to thermal stress weathering. where α, E, ν denote the assumed macroscopic (homogenized) isotropic thermoelastic properties of the rock. z is the rock depth, shown in Fig. 1 and note that ζ is the integral dummy variable. Here, the assumed macroscopic properties are taken to be α " 8.5 ˆ10 ´6 K ´1, E " 74.2 GPa, and ν " 0.28, which are representative of ordinary chondrite materials [START_REF] Delbo | Thermal fatigue as the origin of regolith on small asteroids[END_REF]. Building on the analytic solution for macroscopic thermal stress fields in an isotropic, homogeneous plate, c.f. (Bruno A. Boley, 2012), [START_REF] Wilkerson | Multiscale Mechanics of Failure in Extreme Environments by[END_REF] showed that the corresponding normal component of the macroscopic thermal stress fields in an isotropic, homogeneous sphere is well approximated by Eq. (3). This solution assumes that the temperature field varies spatially with depth and is otherwise constant in the other two spatial dimensions, i.e. T pz, tq. In section 5, it will be shown that non-linear temperature fields, i.e. ∇ 2 T ‰ 0, are required to generate non-zero values of the macroscopic stress Σ. Q in Eq. ( 3) denotes an empirical dimensionless geometrical factor obtained from finite element calculations bridging between thermal stress values in a semi-infinite plate and those of a sphere. Q " 1 for a semi-infinite plate, and varies from Q "1/3 to Q "3/4 for spherical geometries. Informed by finite element calculations, here we approximate this geometrical factor as

Q " max ˆ1 3 , 3 4 D `4δ D `12δ ˙(4)
where δ is the thermal skin depth, which is an approximate measure of the depth at which surface temperature amplitude is attenuated by a factor of 1{e, [START_REF] Titus | thermal diffusivity experiment at the grand falls dune fields[END_REF], i.e.

δ "

d kP πρc p . ( 5
)
Microscopic stress generated by heterogeneities:

σ T M pz; y, tq " 1 2 K ˚∆α∆T pz, tqˆ(6) " r 3 i ´3 ´1 0 ď }y} ď r i 1 2 r 3 i `}y} ´3 ` ´3˘r i ď }y} ď ,
where ∆T pz, tq " 2ˆ`T pz, tq ´T ˘denotes the current diurnal temperature variation at various depths with T denoting the average temperature experienced during a full diurnal cycle. K ˚and ∆α denote parameters that depend on the variation in thermoelastic properties of the heterogeneities in the material microstructure.

K ˚" 12µ m K i K m 3K i K m `4µ m K m `4µ m pK i ´Km qr 3 i ´3
∆α " α i ´αm where µ, K and α are shear modulus, bulk modulus and thermal expansion, respectively, with subscripts i and m referring to inclusion and matrix. Here, the values for these two material parameters are taken as K ˚" 75.3 GPa and ∆α " 1.9 ˆ10 ´6 K ´1, which are representatives of ordinary chondrite materials [START_REF] Medvedev | Determination of the physical properties of stony meteorites with application to the study of processes of their destruction[END_REF][START_REF] Flynn | Physical properties of meteorites and interplanetary dust particles: clues to the properties of the meteors and their parent bodies[END_REF][START_REF] Robert | Temperature, pressure and composition: Structurally analogous variables[END_REF].

Heat diffusion towards the subsurface:

ρc p BT pz, tq Bt " k B 2 T pz, tq Bz 2 , (7)
governs the temperature fields that evolve throughout the day as a function of time t and depth z. The celestial body being analyzed is assumed to be significantly larger than the thermal skin depth, and hence an asymptotic boundary conditions is utilized at large depths, i.e. BT {Bz| z"δ " 0. We consider a surface element (facet) located at the equator of an atmosphereless body whose axis of rotation is perpendicular to its orbit. For such a facet, energy conservation provides the surface boundary:

σ B T 4 pz " 0, tq ´BT pz, tq Bz ˇˇˇz "0 " p1 ´AqS @ r 2 xcos ωty (8)
where x¨y denotes the Macaulay brackets, S @ is the solar constant at 1 AU, A the bolometric Bond's albedo, the infrared emissivity, and σ B the Stefan-Boltzman constant.

Following the arguments of [START_REF] Delbo | Thermal fatigue as the origin of regolith on small asteroids[END_REF], we assume that the temperature field in a rock that is embedded within regolith may be adequately approximated by a one-dimensional thermal transport model. Further details of the assumed radiative boundary conditions utilized are provided in [START_REF] John R Spencer | Systematic biases in radiometric diameter determinations[END_REF]. Here, the thermophysical properties are taken to be representative of those found on S-type asteroids, i.e. thermal diffusivity of k{ρc p " 1.08 mm 2 s ´1, thermal inertia of Γ " 1, 800 J m ´2 s ´0.5 K ´1, and a bolometric albedo of A " 0.1 [START_REF] Opeil | The thermal conductivity of meteorites: New measurements and analysis[END_REF].

Here we assume that k, ρ and c p are depth and temperature independent, which is known to be an approximation. Temperature dependent parameters would slightly increase the values of diurnal temperatures as shown in Fig. 7 of [START_REF] Delbo | Asteroid thermophysical modeling[END_REF].

Together, Eqs. ( 2)-( 8) constitute the necessary set of equations whose numerical solutions provide predictions of the thermal stress of rocks on airless bodies. In the next section, results of this model will be utilized to provide a clearer picture of the key factors governing thermal stress in the inner solar system.

Key factors governing diurnal temperature vari-

ations, rates, gradients, and fields

Scaling diurnal surface temperature variations

We begin with a discussion of the factors influencing the maximum diurnal surface temperature variations (difference between the minimum and the maximum diurnal temperatures), denoted here as ∆ Ts , with various rotation periods P " 2π{ω and orbital distances r. Throughout this paper ˆwill generally denote the maximum value of said variable. According to Eq. ( 6), the magnitude of the maximum microscopic stresses near the surface are directly proportional to ∆ Ts , and hence the diurnal surface temperature variation is expected to be a vital factor governing the thermal stress weathering lifetime of rocks.

Numerical solutions of Eqs. ( 7)-( 8) are shown in Fig. 2 for a subsolar surface element on the equator. The open circles in Fig. 2 demonstrate how the maximum diurnal surface temperature, i.e. ∆ Ts , is affected by the rotation period P on S-type asteroids at r " 1 AU. The rotation period is varied from P " 6 minutes (240 rotations per day) representative of very fast (likely small and monolithic) rotators, e.g. 2000 DO 8 , to P " 29.5 Earth days (708 hours) representative of very slow rotators, e.g. Earth's Moon. ∆ Ts is found to increase monotonically with increasing rotation period, i.e. slow rotators experience larger diurnal surface temperature variations than fast rotators.

Fig. 2: The magnitude of ∆ Ts{ Teq, ∇ Ts, ∇ 2 Ts and d Ts{dt at the surface for a boulder located on the equator at the surface of a body at 1 AU with respect to rotation periods (where ˆgenerally denote the maximum value of said variable in this paper). The markers denote numerical calculations whereas the solid lines show approximations due to Eq. ( 11)-( 12). The labels of the vertical axis are given in the legend of the figure.

The monotonically increasing trend displayed by the open circles in Fig. 2 is to be expected, since in the limit of an infinitely fast rotator, i.e. P Ñ 0, the diurnal temperature variation will approach zero, i.e. ∆ Ts Ñ 0, due to the fact that an infinitely fast rotator will experience nearly isothermal conditions for each latitude. On the other hand, the subsolar surface on a slow rotator is exposed to solar rays for longer periods of time, thereby allowing it to heat to higher temperatures. In the limit of an infinitely slow rotator, the maximum subsolar diurnal temperature variation will approach the maximum subsolar equilibrium model temperature [START_REF] Frank | Fundamentals of heat and mass transfer[END_REF]

, i.e. ∆ Ts Ñ T eq " 4 c 1 ´A σ B S @ r 2 . ( 9
)
For the S-type asteroids modeled here, T eq " 394 K at r " 1 AU. A rotation rate of P " 708 hours is not slow enough for ∆ Ts to converge to T eq , and extrapolation of the results shown in Fig. 2 suggests that this convergence will occur for P Á 10 5 ´10 6 hours.

The thermal parameter is a non-dimensional measure of how effectively a body's surface temperature can keep pace with diurnal changes in the insolation [START_REF] John R Spencer | Systematic biases in radiometric diameter determinations[END_REF], and is expressed as

Θ " Γ ? ωr 3 4 a σ B p1 ´Aq 3 S 3 @ . (10)
Qualitatively, the delay between a change in illumination on a surface element and the corresponding change in temperature, i.e. thermal inertia, can be masked by a very slow rotation. In other words, a surface element on the non-illuminated part of an extremely slow rotator will be cold regardless of the value of thermal inertia because the changes in illumination happen very slowly and one cannot tell whether it is cold because there was no heat conducted towards the subsurface (low thermal inertia) or because there was sufficient time for the energy in lower layers to be radiated into space. Conversely, a facet on an extremely fast rotator will stay warm when it is not illuminated even if the thermal inertia is low. The thermal parameter, Θ, defined as the ratio between the radiative and the rotational time scales, is what determines the effects of heat conduction towards the subsurface on the time evolution of the temperature. Θ " 10 ´1 ´10 2 for most near-Earth asteroids. Isosurface lines of constant Θ as a function of orbital distance and rotational period are provided for reference in Fig. 5. For moderately fast rotators, i.e. P À 10 hours, the maximum diurnal surface temperature variation is found to be roughly proportional to the inverse of the thermal parameter, i.e. ∆ Ts 9Θ ´1, which implies that Θ∆ Ts is approximately constant for all fast rotators. On the other hand, for slow rotators, i.e. ω Ñ 0, ∆ Ts Ñ ζ 8 T eq where ζ 8 " 2{3 defines the fraction of T eq that ∆ Ts is apparently converging to in Fig. 2 for the maximum rotation period of interest, i.e. P " 10 3 hours (i.e. ∆ Ts { Teq is about 2/3 for rotation period of 10 3 hours in Fig. 2). Here, a simple empirical relationship is proposed to smoothly transition between these two regimes, i.e.

∆ Ts

" Θ lim ωÑ8 Θ∆ Ts `1 ζ 8 T eq , (11)
where lim ωÑ8 Θ∆ Ts " 1.12 ˆTeq is a reference scaling factor obtained from extrapolation of the numerical calculations provided in Fig. 2.

Considering its simplicity, the agreement between Eq. (11) (solid black line in Fig. 2) and the numerical calculations (open black circles) is quite remarkable. Even more remarkable is the fact that Eq. ( 11) accurately predicts the effect of varying thermophysical properties as well as orbital distances on maximum diurnal surface temperature variations as shown in Fig. 3. For example, Θ for S-type asteroid is roughly 3 times greater than corresponding case of C-type asteroid with Γ = 640 J m ´2 s ´0.5 K ´1 and A " 0.02 [START_REF] Delbo | Thermal fatigue as the origin of regolith on small asteroids[END_REF], and hence ∆ Ts will be roughly 3 times higher on fast rotating C-type asteroids as compared with fast rotating S-type asteroids at the same spin rate and orbital distance. On the other hand, T eq is roughly equal on both S-type (T eq " 394 K at 1 AU from the Sun) and C-type (T eq " 402 K at 1 AU from the Sun) asteroids, and hence ∆ Ts will be roughly equal for slow rotating S-type and Ctype asteroids at the same spin rate and orbital distance. Likewise, the maximum diurnal surface temperature variation monotonically decreases with orbital distance as r ´2 for fast rotators and r ´1{2 for slow rotators according to Eqs. ( 9)-( 11). That said, it should be noted that in general a fast rotator, from a thermal diffusion standpoint, is defined as a body with Θ " 1, and hence the minimum spin rate necessary to classify a body as a fast rotator depends on the orbital distance as shown in Fig. 4. For example, the minimum spin rate to classify a body as a fast rotator at 1 AU is 8 times smaller at 2 AU according to Eq. ( 10). 11) shown as solid red line, and numerical data shown as blue circles. The numerical data includes predictions at various orbital distances and for 2 sets of thermophysical material properties, i.e. ordinary and carbonaceous chondrites.

According to [START_REF] Molaro | Grain-scale thermoelastic stresses and spatiotemporal temperature gradients on airless bodies, implications for rock breakdown[END_REF], the equatorial diurnal temperature variations on the Moon, i.e. ∆ T " 290 K, may drive thermal weathering. They made it clear that the threshold required to drive crack propagation via thermal stresses is not well constrained in this environment. Nevertheless, if we consider this temperature amplitude as a sufficient driving force on the Moon, it is important to note that same situation is possible at different conditions as well. For example, at 1 AU, a similar equatorial diurnal temperature variation is achieved on very slow rotating S-type asteroids with a period of P Á 640 hours, according to Fig. 2. That said, the same diurnal temperature variation can be achieved on much faster rotating S-type asteroids whose orbits bring them closer to the Sun, e.g. P Á 25 hours at 0.5 AU and P Á 1.4 hours at 0.1 AU.

Scaling time-rate-of-change in surface temperatures

from diurnal temperature variations

In addition to diurnal temperature variation, the maximum time-rate-of-change in surface temperatures, i.e. |d Ts {dt|, has long been considered (albeit with some confusion) an even more important factor governing the efficacy of thermal stress weathering [START_REF] Richter | Thermal expansion behavior of igneous rocks[END_REF]. Numerical solutions of the maximum time-rate-of-change in surface temperatures as a function of rotation period (at r " 1 AU) are shown as red stars in Fig. 2 for S-type asteroids. Here, |d Ts {dt| decreases monotonically with increasing rotation period P , i.e. fast rotators have higher surface temperature rates than do slower rotators. This trend is opposite from that observed for ∆ T , which scales roughly as P 1{2 on fast rotators, where as |d Ts {dt| scales roughly as P ´1{2 on fast rotators. An overly simplified, yet illustrative, explanation for this particular scaling is that the average surface temperature time-rate-of-change is proportional to the ratio of the diurnal temperature variation and the period over which the temperature changes from maximum to minimum, i.e. |dT s {dt| avg 9∆ T { P 2 9P ´1{2 . Hence, provided that |d Ts {dt|9|dT s {dt| avg then such a scaling relationship is expected. On the other hand, with respect to orbital distance r, |d Ts {dt| scales very similarly to that of ∆ T as shown in Fig. 4. Therefore, we make use of these arguments to provide a simple relationship for the dependence of the maximum time-rate-of-change in surface temperature on rotation period and orbital distance, i.e.

ˇˇd Ts {dt ˇˇ« β ˇˇd Ts {dt ˇˇa vg « 2β ∆ Ts P , (12)
where ∆ T is approximated via Eq. ( 11) with β « 2.5 being a proportionality constant between |d Ts {dt| and |dT s {dt| avg .

Remarkable agreement between numerical calculations of |d Ts {dt| (shown as pink stars) and the simple approximation given by Eq. ( 12) (pink lines) is shown in Fig. 2.

While it has received warranted criticism, c.f. Boelhouwers and Jonsson (2013); [START_REF] Viles | Simulating weathering of basalt on Mars and Earth by thermal cycling[END_REF] and Eppes et al. (2016), in much of the literature |d Ts {dt| Á 2 K/min is assumed to be a threshold criterion that must be achieved in order for thermal stress weathering to be an effective process on bodies with thermal diffusivities similar to rocks [START_REF] Hall | The role of thermal stress fatigue in the breakdown of rock in cold regions[END_REF]. For the S-type asteroids modeled here, this criterion is met for very fast rotating near-Earth asteroids, e.g. P À 1 hour at r " 1 AU, or moderately fast rotating asteroids with very small perihelion, e.g. P À 10 hours at r " 0.1 AU.

Consider now that either |d Ts {dt| Á 2 K/min or ∆ T Á 100 K could plausibly be taken as overly simplified criterion for the initiation of thermal stress weathering. Here we compute spin and orbital distance conditions that met either criteria. In this simple case, thermal stress weathering would be anticipated on either fast or slow rotators as shown in Fig. 5 regions I, II and III. For near-Earth S-type asteroids, thermal stress weathering would be expected to be ineffective on asteroids with intermediate rotation rates, e.g. 1 hr À P À 100 hr , region IV . This band-gap will increase with increasing perihelion distance, e.g. 0.1 hr À P À 10 3 hr for S-type asteroids in the Main-belt (r " 2.5 AU). Likewise, the width of this band-gap will decrease with decreasing perihelion, and completely vanish at sufficiently small perihelions, e.g. r À 0.5 AU for S-type asteroids, as shown in Fig. 5. A vanishing band-gap implies that at such small orbital distances, thermal stress weathering would be anticipated to be effective on all S-type asteroids, regardless of spin rate.

Noting correlations in temporal and spatial gradients of surface temperatures

While Fig. 5 is helpful to elucidate some of the key governing factors affecting thermal stress weathering, the criterion utilized are overly simplistic, and will be shown to be insufficient in subsequent sections of this paper. In particular, we note that the temporal gradients of surface temperature, i.e. |d Ts {dt|, do not directly enter into the calculation of macroscopic stress, i.e. Eq. (3). Rather, the macroscopic stress distribution through the interior of the rock, i.e. Σpz, tq, is a function of the degree of non-linearity in the spatial distribution of temperature T pz, tq. The magnitude Fig. 5: A map for S-type asteroids indicating the spin rates and orbital distances necessary to achieve classical criteria for thermal stress weathering, i.e. ∆Ts ě 100 K and dTs{dt ě 2 K/min, according to Eq. ( 11) and Eq. ( 12). Isosurface lines indicate corresponding values of the thermal parameter Θ on S-type asteroids.

of the second-order spatial gradient of surface temperature, i.e. |∇ 2 Ts |, is one (albeit insufficient) measure of the degree of non-linearity in T pz, tq.

Numerical calculations of the dependence of |∇ 2 Ts | on rotation period P are shown (as blue triangles) in Fig. 2, where the dependence on P is identical to that of |d Ts {dt|. Such a scaling should be anticipated, since these quantities are proportional according to the heat equation, i.e. |d Ts {dt| " pk{ρc p q|∇ 2 Ts |. Therefore, |∇ 2 Ts | may be adequately approximated as the quotient of Eq. ( 12) and thermal diffusivity, which is shown (as a blue line) to agree with numerical calculations (blue triangles) in Fig. 2.

Additionally, the scaling of |∇ 2 Ts | with orbital distance is shown in Fig. 4. Note that the assumed thermal stress weathering criterion of |d Ts {dt| Á 2 K/min is equivalent to |∇ 2

Ts | Á 5.5 K/cm 2 for S-type asteroids. Obviously, Ts | will be shown to be an insufficient measure of how macroscopic surface stresses scale with rotation period and orbital distance. For completeness, Fig. 2 and Fig. 4 provide numerical calculations (red diamonds) of |∇ Ts | as a function of rotation period and orbital distance, respectively. In subsequent sections, it will be shown that (as should be expected) the first spatial gradient of surface temperature, |∇ Ts |, is not correlated with thermal stresses as noted by [START_REF] Molaro | Grain-scale thermoelastic stresses and spatiotemporal temperature gradients on airless bodies, implications for rock breakdown[END_REF].

|d

Effect of diurnal temperature variations and gradients on thermal stress fields

Here we turn our attention to elucidating the manner by which diurnal temperature variations and spatial gradients govern microscopic and macroscopic thermal stresses. Section 5.1 provides a discussion of microscopic (grain-scale) thermal stresses and their correlation with the diurnal temperature variation ∆T . Our analytic analysis largely agrees with the overall trends of and attempts to generalize the two-dimensional finite element calculations of [START_REF] Molaro | Grain-scale thermoelastic stresses and spatiotemporal temperature gradients on airless bodies, implications for rock breakdown[END_REF], although the magnitude of our calculated microscopic thermal stresses are roughly an order of magnitude lower than those of [START_REF] Molaro | Grain-scale thermoelastic stresses and spatiotemporal temperature gradients on airless bodies, implications for rock breakdown[END_REF] which made use of rigid boundary conditions. Section 5.2 provides a discussion of macroscopic (rock-scale) thermal stresses and their complex dependence on rock size, thermal skin depth, and second-order (and higher) spatial gradients of the temperature field. Our analytic analysis largely agrees (both in overall trends and magnitudes) with the full three-dimensional finite element calculations of [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF], and in subsection 5.3 builds upon and attempts to generalize their results through simple scaling laws that can be easily applied to a wider range of spin rates and orbital distances as well as to bodies with different thermophysical and thermoelastic material properties. Lastly, subsection 5.4 generalizes and provides a physical interpretation of the non-monotonic dependence of macroscopic thermal stresses on rotational period first reported by [START_REF] Ravaji | The effect of rotation period on thermal stress weathering[END_REF] Note that throughout the following sections, variation indicates the difference between the maximum and minimum values of a variable experienced throughout a full thermal cycle at a given spatial location. The term gradient corresponds to a local spatial derivative of a particular field quantity.

Proportionality between diurnal temperature variations and microscopic thermal stresses

Microscopic thermal stresses (σ T M ) develop due to spatial fluctuations in the local thermal expansion coefficients, i.e. ∇α ‰ 0, throughout the heterogeneous structure of rocks. For bodies heated primarily by solar radiation, the maximum microscopic thermal hoop stress, i.e. σT M s , will occur near the surface where temperature variations are largest, and may be computed by maximizing Eq. ( 6), i.e.

σT M s " K ˚∆α∆ Ts max ˆ1 ´vi , 1 `vi 2 ˙, (13)
where v i " r 3 i { 3 is the volume fraction of chondrules. Due to the proportionality between diurnal surface temperature variations and maximum microscopic thermal stress, i.e. σT M s 9∆ Ts , spin rate, orbital distance, and thermophysical material properties have identical affects on both maximum microscopic thermal stresses and diurnal surface temperature variations. For example, the effect of rotational period P on σT M s will follow the same monotonically increasing (non-linear) trend as the open black circles in Fig. 2, as demonstrated in [START_REF] Ravaji | The effect of rotation period on thermal stress weathering[END_REF]. Likewise, the monotonically decreasing (non-linear) dependence of orbital distance r on σT M s will follow the trend of the open black circles in Fig. 4. The combined effect of rotational period and orbital distance on the maximum microscopic thermal stresses are shown as a contour map in Fig. 6. The contour map is generated via combination of Eq. ( 13) and Eq. ( 11). As expected, microscopic thermal stresses are highest on slow rotating bodies, e.g. Moon, that pass close to the Sun. Note that a diurnal temperature variation of ∆ Ts " 100 K (blue line in Fig. 5) generates a maximum microscopic thermal stress of σT M s " 11 MPa (yellow region in Fig. 6) on S-type asteroids. Additionally, it is evident in comparing Fig. 5 and Fig. 6 that constant values of Θ produce constant values of σT M s .

Recall that for fast rotating asteroids the diurnal surface temperature variation is inversely proportional to the thermal parameter, i.e. ∆ Ts 9Θ ´1, and hence σT M s 9Θ ´1. It so happens that K ˚∆αp1 `vi q and T eq are roughly equal on both S-type and C-type asteroids, where as Θ is roughly 3 times higher on S-type asteroids than C-type asteroids. As such, it follows that the maximum microscopic thermal stresses on slow rotating C-type and S-type asteroids will be roughly equal, where as σT M s will be roughly 3 times higher on fast rotating C-type asteroids.

Microscopic stresses were also analyzed by [START_REF] Molaro | Grain-scale thermoelastic stresses and spatiotemporal temperature gradients on airless bodies, implications for rock breakdown[END_REF], who reported thermal mismatch stresses on the order of 100-200 MPa for rocks on the Moon. This discrepancy is primarily associated with a different choice of boundary conditions. In their case the material is assumed to be rigidly confined in both lateral directions. Such boundary condition causes that even a homogeneous body will develop (in-plane equibiaxial) thermal stresses under uniform temperature field, i.e. σ " ´Eα∆T . However, in [START_REF] Delbo | Thermal fatigue as the origin of regolith on small asteroids[END_REF] traction-free boundary conditions are utilized to represent an unconfined rock resting on a near frictionless surface. Under these boundary conditions, a homogeneous body, i.e. ∆α " 0, subject to uniform temperature fields, i.e. ∇T " 0, will not develop micro scale stresses. These later boundary conditions more closely mimic the reality of a rock undergoing thermal expansion while resting upon a bed of regolith, since the boulder is free to expand and contract on the surface. Even when the boulder is embedded within regolith, it is doubtful that the regolith can impose any significant constraints, since the rock's stiffness is orders of magnitude larger than that of the regolith. To make our results comparable, the effect of this rigid boundary condition may be approximately removed from [START_REF] Molaro | Grain-scale thermoelastic stresses and spatiotemporal temperature gradients on airless bodies, implications for rock breakdown[END_REF] reported results through the principle of superposition, i.e. σ " σ `σT M , where σ is the average stress in the domain. The reported values of σ and σ lead to σ T M " 36 MPa for P " 700 hr at 1 AU. The primary difference between this value and ours is associated with the different materials assumed in each study, K ˚" 134 GPa and ∆α " 0.4 ˆ10 ´5 K ´1 in [START_REF] Molaro | Grain-scale thermoelastic stresses and spatiotemporal temperature gradients on airless bodies, implications for rock breakdown[END_REF] versus K ˚" 75 GPa and ∆α " 0.19 ˆ10 ´5 K ´1, which explains the roughly factor of 4 discrepancy between the modified values from σ T M " 36 MPa in [START_REF] Molaro | Grain-scale thermoelastic stresses and spatiotemporal temperature gradients on airless bodies, implications for rock breakdown[END_REF] and our calculations (σ T M " 9 MPa) for P " 700 hr at 1 AU. Given the various assumptions, the agreement between our scaling relations for microscopic thermal stresses and the modified finite element results of [START_REF] Molaro | Grain-scale thermoelastic stresses and spatiotemporal temperature gradients on airless bodies, implications for rock breakdown[END_REF] is quite remarkable. This agreement provides a degree of confidence in the scaling relationships, which are far more general and accessible for the broader community.

The complex relationship between macroscopic stresses and spatiotemporal temperature fields

We now turn our attention to the spatiotemporal macroscopic (in-plane equibiaxial) thermal stresses, i.e. Σpz, tq, and its complex relationship with the spatiotemporal temperature fields, i.e. T pz, tq. Here the spatiotemporal temperature and macroscopic stress fields are computed via Eq. ( 7) and Eq. (3), respectively. Examples of T pz, tq and Σpz, tq as well as dT pz, tq{dt are shown in Fig. 7 for a 10 cm boulder on a fast rotating (P " 0.6 hr) and a slow rotating (P " 103 hr) S-type asteroid at 1 AU from the Sun. Our results are qualitatively similar to that reported by [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF] despite our different modeling approach and choice of thermomechanical properties, rotation periods, and boulder sizes. Comparing Fig. 7(a) and (b), it is evident that the peak tensile macroscopic stress in the faster rotator (P " 0.6 hr) is Σ " 0.7 MPa, which is roughly 2.5 times larger than the peak macroscopic stress in the slow rotator (P " 103 hr) of Σ " 0.3 MPa. This trend is opposite of what was observed for microscopic stresses, i.e. σT M s " 0.9 MPa and σT M s " 7 MPa for P " 0.6 hr and P " 103, respectively, at 1 AU, c.f. Fig. 6. As discussed earlier, the peak microscopic stresses are proportional to temperature variation, ∆T s , which increases with increasing rotation period, as shown in Fig. 2. On the other hand, the peak macroscopic stress, Σ, is more closely correlated with second gradient of temperature, ∇ 2

Ts , which was shown in Fig. 2 to decrease with increasing rotation period.

This point is further clarified by comparing the macroscopic surface stress reported in Fig. 7(a) and (b) with corresponding surface temperature in Fig. 7(c) and (d), for which it is evident that the peaks of temperature history and macroscopic thermal stress are out of sync. In the case of the slow rotator (P " 103 hr), thermal stress is negligible, Σ " 0, when the temperature is at its maximum and minimum values. This is due to the fact that the macroscopic stress is a function of temperature gradients and not the temperature magnitude. The maximum and minimum temperatures occur when temperature rate at depth z is zero, dT {dtpzq " 0, c.f. Fig. 7(d) and (f), which for a slow rotator implies that the second-order thermal gradients at depth z are small, i.e. ∇ 2 T pzq " 0 via the heat conduction equation, and hence the macroscopic stresses are also small. On the other hand, the peak macroscopic stress at various depths, i.e. Σpzq, is well correlated with the the peaks of dT {dtpzq as observed in comparing Fig. 7 For slow rotators, a remarkable agreement between the shape of the surface temperature time-rate-of-change profile, i.e. solid lines in Fig. 7(f), and the macroscopic stress, i.e. Fig. 7(b), is observed. This implies that spatiotemporal macroscopic stress is proportional to rate of temperature change at surface, i.e. Σpz, tq9dT s ptq{dt ˆf pzq, where f pzq determines the spatial dependence of the macroscopic stress. It is interesting to note that for z " D{5 in Fig. 7(b) that the macroscopic stress are nearly zero throughout the day. Although not shown, the macroscopic stresses at z " 4D{5 are also roughly zero throughout the day. Furthermore, the macroscopic stress will tend to reverse sign at these locations, e.g. a tensile stress for z À D{5 and z Á 4D{5, but a compressive stress for D{5 À z À 4D{5, and vicea-versa. This is somewhat analogous to mechanical bending with z " D{5 and z " 4D{5 acting as the so-called neutral axes. The sign of the thermal stress at the surface can generally be determined via the following procedure. For slow rotators, the surface stress will be tensile from sunset to sunrise (cooling phase), and compressive from sunrise to sunset (heating phase). The interior (D{5 À z À 4D{5) thermal stress will generally follow the opposite trend, e.g. compressive from sunset to sunrise, and tensile from sunrise to sunset, c.f. Fig. 7(b), (d), and (f). For fast rotators, the surface macroscopic stress is generally compressive when the surface temperature is higher than the temperature at some depth, e.g. z " D{5. Likewise, when the surface temperature is cooler than at z " D{5, the surface stress tends to be in tension. For example, in the beginning of the solar day, on the surface, temperature is higher than that at D{5 and therefore the corresponding thermal stress is negative (compressive) in Fig. 7(a). If we roughly assume z " D{5 " 2 cm is one of the neutral axes of the 10 cm rock in Fig. 7 (shown as yellow curve in the plots). In other words, the surface stress changes sign when the surface temperature crosses over the yellow curve in Fig. 7(c). Similar trends are also reported by Eppes et al. (2016) where observed cracking, as measured by acoustic emissions, occurs at the same times of the day that our modeling predicts.

In comparing Fig. 7(c) and (d), we can see one of the fundamental differences between the fast rotators and slow rotators when it comes to macroscopic thermal stresses. In this case of the slow rotator shown Fig. 7(d) and (f), the temperature history and its time derivative do not vary too much with respect to location. Hence, the corresponding maximum, minimum, and zero thermal stresses occur almost at nearly the same time, regardless of depth, as shown in Fig. 7(b). However, for the fast rotator shown in Fig. 7(c), there is a significant phase lag (roughly z{δ) between temperature on the surface and temperature at other depths, which complicates matters. In this case, the temperature history and its time derivative vary significantly with respect to location and thus, the peak of the corresponding thermal stress shifts as well. Furthermore, due to this phase lag, the shape of the macroscopic thermal stress profiles are completely different than the shape of the temperature rate profiles, c.f. Fig. 7(a) and (e), whereas the shapes of the stress and temperature rate profiles are nearly identical on slow rotators, c.f. Fig. 7(b) and (f).

To better understand the relationship between temperature fields and stress fields, it is illustrative to consider an idealized temperature field whose spatial dependence can be described by an n-th order polynomial, i.e. T pzq " ř n j"0 c j z j , for which the integrals in Eq. ( 3) may be obtained analytically, i.e.

Σpzq " QαE p1 ´νq

n ÿ j"2 " 6jc j D j´1 z pj `1qpj `2q ´cj z j (14) ´2pj ´1qc j D j pj `1qpj `2q * .
For a quadratic temperature field, i.e. T pzq " the stress field further simplifies, i.e.

c 0 `c1 z `c2 z 2 , (a) (b) (c) (d) (e) (f)
Σpzq " Σ s ˆ1 ´6 z D `6 z 2 D 2 ˙, (15)
where Σ s " ´QαED 2 c 2 {6p1 ´νq. The spatial function in the parenthesis in Eq. ( 15) is symmetric for the domain 0 ď z ď D, and has a maximum value of 1 at z " 0 and z " D, a minimum value of ´1{2 at z " D{2, and a zero-value at z " D{5 and z " 4D{5 (i.e. the bending neutral axes), which is in fairly good agreement with the spatial dependence of the macroscopic stress shown in Fig. 7(b). In general, an n-th order spatial dependence of the temperature field will lead to an n-th order spatial dependence of the macroscopic stress field (with n bending neutral axes).

It is interesting to note that c 0 and c 1 are exactly canceled out when analytically carrying out the integration in Eq. ( 14). Physically, c 0 represents the uniform component of the temperature field and c 1 is the linearly varying component of the temperature field, neither of which apparently affect macroscopic thermal stresses. Alternatively, one can express these coefficients as the surface temperature and the n-th order spatial derivatives of the surface temperature, i.e. c 0 " T s and c n " ∇ n T s . Using this expression of the coefficients, the macroscopic thermal stresses at the surface of the boulder may be expressed as

Σ s " ´QαE p1 ´νq ˆD2 6 ∇ 2 T s `D3 5 ∇ 3 T s `... (16)
`2pn ´1qD n pn `1qpn `2q

∇ n T s ˙.
Thus, the macroscopic thermal stresses (of an unconfined body) are independent of both the surface temperature and its first spatial derivative. Moreover, it is clear from Eq. ( 16) that negative (positive) thermal gradients are associated with tensile (compressive) surface stresses.

Clearly the relationship between macroscopic thermal stresses and the spatial gradients of the temperature field is rather complex. Deeper insights will provided in section 5.3 for a few simpler, yet illustrative, special cases.

Simple scaling relations for macroscopic stresses

Given that second-order and higher-order spatial derivatives of surface temperature will govern the macroscopic thermal stresses, it is convenient to obtain an approximate scaling relation for the higher-order spatial gradients. A rough approximation of the higher-order spatial gradients of surface temperature may be obtained by considering the simplest case in which the diurnal surface temperature cycle follows simple harmonic oscillation, i.e. T s ptq " Ts `1 2 ∆ Ts cospωtq, where Ts is the mean surface temperature. This simple temporal variation of the surface temperature gives rise to a spatiotemporal temperature field below the surface that can be obtained analytically from Eq. ( 7), c.f. [START_REF] Anderson | Near-Surface Thermal Profiles in Alpine Bedrock: Implications for the Frost Weathering of Rock[END_REF], i.e.

T pz, tq " Ts `1 2 ∆ Ts e ´z δ cos ´ωt ´z δ

¯. ( 17
)
Following Eq. ( 17) the maximum magnitude of the secondorder and higher-order thermal gradients are approximated as ˇˇ∇ n Ts ˇˇ" ?

2 n´2 ∆ Ts δ n . ( 18
)
Through inspection of Eq. ( 16) and Eq. ( 18), it is evident that each term in Eq. ( 16) corresponding with an nth degree temperature gradient will roughly scale with the product of diurnal surface temperature variation and the n-th power of the ratio of rock diameter to thermal skin depth, i.e. pD{δq n ∆ Ts . As such, second order derivatives of the temperature field will tend to govern the magnitude of macroscopic thermal stresses on rocks that are small in relation to the thermal skin depth, i.e. Ts ˇˇ" pρc p {κq ˇˇd Ts {dt ˇˇ" ∆ Ts {δ 2 denoting the maximum magnitude of the second-order spatial gradient of the surface temperature. The importance of second order thermal gradients (and the insignificance of first order thermal gradients) is emphasized in Eq. ( 19). The effect of material properties on macroscopic surface stresses is evident. For example, consider two similar sized (small) rocks that are undergoing the same time-rate-of-change in surface temperatures dT {dt on a slow rotating body. In this case, the macroscopic thermal stresses at the surface will vary linearly with thermal expansion coefficient α and effective stiffness E{p1 ´νq, while varying inverse linearly with thermal diffusivity κ{ρc p . On the other hand, if one considers two small rocks of different composition on a sufficiently slow rotating body (such that ∆ Ts Ñ T eq ) then the macroscopic thermal stresses at the surface will scale roughly with 4 a p1 ´Aq{ (in addition to the previously mentioned scalings). These scaling dependencies are fairly consistent with the rock lithology sensitivity study (with D " 1 m and δ " 0.8 m) reported in [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF]. [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF] carried out fully three-dimensional finite element calculations of macroscopic thermal stress development in various sized rocks on the Moon (δ " 0.8 m), and reported values of the maximum surface principal stresses of "2 MPa, "5 MPa, and "9 MPa for D " 0.3 m, D " 0.5 m, D " 0.7 m. Remarkably, these stress values are fairly consistent with predictions (respectively "1.5 MPa, "4 MPa, and "8 MPa for the same set of material properties) of Eq. ( 19) despite its simplicity and lack of fitting parameters. These stress values are generally smaller than fracture toughness of typical rocks [START_REF] Gui | Modelling the dynamic failure of brittle rocks using a hybrid continuum-discrete element method with a mixed-mode cohesive fracture model[END_REF][START_REF] Eppes | Mechanical weathering and rock erosion by climate-dependent subcritical cracking[END_REF][START_REF] Asadi | Investigation of Hydrodynamically Dominated Membrane Rupture, Using Smoothed Particle HydrodynamicsFinite Element Method[END_REF][START_REF] Mir | The efficiency of thermal fatigue in regolith generation on small airless bodies[END_REF] and would to slow steady fatigue process.

Σs

The slight under-prediction of Eq. ( 19) in comparison with the finite element calculations of [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF] is likely due to (i) the fact that Eq. ( 17) under-predicts the timerate-of-change in surface temperature near sunrise and sunset; (ii) the three-dimensional nature of the finite element calculations; and (iii) that [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF] reports the maximum principal stresses whereas we report the in-plane normal stress components. The maximum principal stress and the maximum in-plane normal stress reported here are equivalent in the absence of shear stresses; however, the maximum principal stress will be larger if shear stresses are nonnegligible. 2017) for a boulder on the Moon surface with analytical approaches provided here via Eq. A.3 and Eq. ( 21) for the maximum stress at the surface, Σs, and Eq. ( 22) for half rock diameter, ΣD{2 , with material properties adopted from [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF]. Although stress magnitudes are different (most likely because of the difference between 1D and 3D conditions), trend of the plots are in fairly good agreement.

In addition, [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF] carried out threedimensional finite element calculations on larger rocks, and observed that the magnitude of surface stresses increases sublinearly with rock size for D ą 10 m as shown in Fig. 8. Since D ą δ " 0.8 m for these finite element calculations, Eq. ( 19) is not applicable in this regime. Further insights into size effects in this regime may be gleaned by considering temperature fields for which the phase lag in Eq. ( 17) is neglected, i.e. T pzq " Ts `1 2 ∆ Ts exp ´z{δ. The macroscopic stress field corresponding to this temperature field may be obtained analytically through integration of Eq. (3), i.e.

Σpzq "

QαE∆ Ts 6p1 ´νq

ˇˇˇe ´z δ ´δ D ´1 ´e´D δ ¯(20) `3 δ D ´2 z D ´1¯" 2 ´ˆ2 δ D `1˙´1 ´e´D δ ¯*ˇˇˇˇ.
The significance of the ratio of rock diameter to thermal skin depth is again apparent in Eq. ( 20). Additionally, the spatial variation of the maximum macroscopic thermal stress has both an exponential decay term, i.e. exp ´z{δ, as well as a linearly varying term p2z{D´1q, which reverses sign at the center of the rock and whose magnitude is greatest near the top (z " 0) and bottom (z " D) of said rock. It is worth noting that such a temperature field will induce an exponential decay in ∆ T pzq and hence the microscopic thermal stresses as well, i.e. σT M pzq9 exp ´z{δ via Eq. ( 7). Similarly, due to the exponential decay term, one can anticipate the magnitude of macroscopic thermal stresses to be greatest near the top of the rock as shown by the solid blue lines in Fig. 7. However, the linearly varying term in Eq. ( 20) can cause the magnitude of macroscopic thermal stresses to increase with increasing depth as shown by the z " 3 cm and z " 5 cm cases in Fig. 7 for D " 10 cm and P " 103 hours. Furthermore, since the linearly varying term reverses sign with depth, this causes part of the rock to be in tension and part to be in compression (as shown in Fig. 7), which may be thought of as being analogous to bending-induced normal stresses in beams [START_REF] Hadi | Size-Dependent Vibration Analysis of FG Microbeams in Thermal Environment Based on Modified Couple Stress Theory[END_REF].

Following Eq. ( 20) the maximum macroscopic thermal stresses at the surface pz " 0q may be conveniently expressed as

Σs " Σ8 s Q Q 8 " 1 ´6 δ D `δ D ˆ6 δ D `2˙´1 ´e´D δ ¯* .
(21) with Q 8 " 3{4 being the geometric factor for D " δ and the macroscopic surface stress for D " δ denoted as Σ8

s " βαE∆ Ts {16p1 ´νq with β being a non-dimensional fitting constant artificially introduced here to bring Σs into closer agreement with higher fidelity calculations (see Appendix A). In the absence of β, Σ8

s " 9.6 MPa for a boulder on the Moon (∆ T " 235 K), whereas [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF] reported a value of "42 MPa, which is a slightly more severe underprediction than our simple model for the special case of D ă δ, i.e. Eq. ( 19). The need for β is associated with the fact that the phase lag term, i.e. cospωt´z{δq in Eq. ( 17), induces larger thermal gradients (and larger time-rate-of-changes in surface temperatures, c.f. Eq. ( 12)). For example, the n´th degree gradients in surface temperature are a factor of ? 2 n (see Eq. ( 18)) larger with the phase lag term than without it. β " 5 brings our simple approximation of macroscopic thermal stresses into better agreement with that of [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF] for D " δ. That said, a value of β " 24 will bring the second-order Maclaurin series (with D{δ Ñ 0) of Eq. ( 21) into agreement with Eq. ( 19). As such, here we choose a simple exponential decay function to transition between these two limits, i.e. β " 24 ´19 expp´δ{Dq.

While this simple analysis fails to fully predict the precise value of Σ8 s , it is rather remarkable that the general size dependence of Σs is captured quite well by Eq. ( 21), as demonstrated by the good agreement with the finite element calculations of [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF] shown in Fig. 8. Furthermore, Fig. 8 includes a higher-fidelity analytic approximation, i.e. Eq. A.3 that directly accounts for the phase lag in Eq. ( 17), of the size-dependence of Σs .

As oppose to increasing stress on the surface, Σs , with respect to the rock diameter, interior thermal stress is not always increasing. Initially, the interior stress is observed to increase with rock size for D{δ À 7, followed by a subsequent decrease in interior stress with increasing rock size for D{δ Á 7. [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF] observed similar behav-ior for the interior maximum stress. Substituting z " D{2 in the Eq. ( 20) results in Eq. ( 22). Stress at the half rock size is shown in Fig. 8 (black line) which captures this nonmonotonic size dependence well.

Σpz " D 2 q " QβαE∆ Ts 6p1 ´νq ˇˇˇ" e ´D 2δ ´δ D ´1 ´e´D δ ¯*ˇˇˇˇ. (22)
A comparison between different approximate solutions provided here, Fourier series approximate solutions provided in Appendix A, and results derived from the numerical approach for the surface stress in a 10 cm rock on an S-type asteroid is shown in Fig. 9 as a function of rotation period. Eq. ( 19) agrees well with the numerical calculations, and (by construction) with Eq. ( 20) for δ Á D. As expected, the 2nd order Fourier series agrees best with the numerical calculations and captures the non-monotonic dependence on rotation period that will be discussed further in the next subsection.

Fig. 9: Comparison of numerical calculations of the maximum macroscopic surface stress Σs for a 10 cm rock on an S-type asteroid at 1 AU with four analytic approximations: a 1st order Fourier series approximation, i.e. Eq. A.3 with n " 1, a 2nd order Fourier series approximation, i.e. Eq. A.3 with n " 2, a 1st order Fourier series approximation with the phase lag neglected, i.e. Eq. ( 20), and a special solution for sufficiently slow rotators such that δ Á D, i.e. Eq. ( 19) . Despite their simplicity, the approximate solutions obtained here for macroscopic thermal stresses at a rock's surface, e.g. Eq. ( 19) (for D ă δ) and Eq. ( 20), are remarkably accurate and provide some important general insights. For example, Eq. ( 19) demonstrates that macroscopic thermal stresses are primarily governed by second-order spatial derivatives of the surface temperature when D À δ. On the other hand, when D " δ, the macroscopic thermal stresses are largely governed by the ratio of thermal skin depth to rock diameter, i.e. δ{D. The former condition holds for sufficiently small rocks on sufficiently slow rotating bodies and the latter condition holds for sufficiently large rocks on sufficiently fast rotating bodies. This holds for both interior and exterior stresses.

Although dependency of thermal stress on the rock size in Fig. 8 may be perceived as a size-effect, from our analysis here, it is clear that this is not a true size-effect in the traditional sense, c.f. size effects in nano-indentation, nanowires, or nanoporous materials, for example. Instead, it is evident here that the ratio δ{D governs the fraction of rock that undergoes a significant diurnal temperature variation, with the remaining (deeper) part of the rock remaining in a nearly isothermal condition. The nearly isothermal portion of the rock (z Á δ) does not experience thermal expansion, while the remaining portion (z À δ) experiences thermal expansion and contraction at sunrise and sunset, respectively. The thermal strain incompatibility between these two regions induces a macroscopic thermal stress that is somewhat analogous to the mechanics underlying microscopic thermal stresses. In the analogy, δ{D is somewhat similar to the volume fraction of chondrules v i in Eq. ( 13). Moreover, the functional form of our model for σT M s , i.e. Eq. ( 13), is remarkably similar to the functional form of Σ8 s with both involving the product of a similar set of thermoelastic material properties, the diurnal surface temperature variation ∆ Ts , and a factor related to the fraction of distinct partitioned regions, i.e. δ{D and v i . In the next section, we will utilize these scaling relations to better understand how macroscopic thermal stresses vary throughout the solar system.

Non-monotonic dependence on rotation period

In section 5.1, the microscopic thermal stresses were shown to increase monotonically with increasing rotation period. On the other hand, [START_REF] Ravaji | The effect of rotation period on thermal stress weathering[END_REF] demonstrated that macroscopic thermal stresses exhibit a non-monotonic dependence on rotation rate. For a D " 10 cm boulder on an S-type asteroid at 1 AU, macroscopic thermal stresses at the surface, i.e. Σs , were shown to increase with increasing rotation period for P À 6 hr and to decrease with increasing period for P Á 6 hr, as shown here by the black markers in Fig. 10. This observation is somewhat peculiar considering that both diurnal surface temperature variations, i.e. ∆ Ts , as well as spatiotemporal gradients in surface temperatures, i.e. d Ts {dt and ∇ 2

Ts , all vary monotonically with spin rate, as shown in Fig. 2. Here we aim to provide fundamental insights into this non-monotonic dependence of rotation period on macroscopic stresses.

The markers in Fig. 10 show numerical calculations of the non-monotonic dependence of Σs on rotation period P for three rock sizes: D " 2 cm (blue diamonds), D " 10 cm (black squares), and D " 30 cm (red circles). In each case, Fig. 10 first increases with increasing rotation period, followed by a subsequent decrease with further increasing rotation period. For each rock size, there is a critical intermediate rotation period P cr for which the macroscopic thermal stresses at the surface are greatest. According to our numerical calculations in Fig. 10, P cr "25 hours, "4 hours, and "0.1 hours for D " 30 cm, 10 cm, 2 cm, respectively. Assuming thermal stress weathering is dominated by macroscopic thermal stresses, one would expect thermal stress weathering to be fastest for rotation periods near these size-dependent critical spin rates.

The question that naturally arises is what factors govern this critical spin rate? Following our discussion at the end of section 5.3, one could postulate that as the rock size becomes larger, macroscopic thermal stresses tend to be greatest when the thermal wave penetrates roughly halfway through the rock, i.e. δ cr " D{2. In this case, the bottom half of the rock experiencing a near isothermal condition would provide a resistance to the desired thermal expansion/contraction of the other half of the rock. Since the two regions are partitioned equally, each region would have roughly the same compliance. As such, neither region would be able to dominate the equilibrium configuration, thus necessitating the maximum possible stresses to overcome the incompatibility of the similarly compliant regions. For the S-type asteroids considered here, the critical thermal skin depths corresponding to the calculated critical spin rates are δ cr « 33 cm, 28 cm, 22 cm, 18 cm, 14 cm, 7 cm, 3.5 cm, and 1 cm for D " 60 cm, 50 cm, 40 cm, 30 cm, 20 cm, 10 cm, 5 cm, and 2 cm, respectively, which are in remarkable agreement with our δ cr " D{2 postulation.

Since the thermal skin depth increases with the square root of the rotation period, i.e. δ9 ? P , the thermal skin depth will be larger than the rock diameter, i.e. δ ą D, for spin rates that are a factor of 4 slower than the critical spin rate, i.e. P Á 4 ˆPcr . Recall that Eq. ( 19), i.e. Σs " αED 2 |∇ 2 Ts |{18p1 ´νq, is an applicable scaling relation for δ ą D. As such, the numerical calculations of Σs on slow rotators with periods of P Á 4 ˆPcr (all markers in Fig. 10 that are a bit to the right of the peak stress for said rock size) should roughly scale with Eq. ( 19). In this regime, we should anticipate macroscopic thermal surface stresses to scale with the square of rock size, i.e. Σs 9D 2 , which the size-dependence of numerical calculations shown in Fig. 10 seem to cooperate. Moreover, the macroscopic thermal stresses in this regime are anticipated to roughly scale with |∇ 2

Ts |, which itself was shown in Fig. 2 to scale roughly with P ´1{2 for sufficiently fast rotators and 1{P for sufficiently slow rotators (such that ∆ Ts Ñ T eq q. The scaling of |∇ 2

Ts | with rotation period alone governs the spin rate dependence of surface thermal stresses for sufficiently slow rotators with P Á 4 ˆPcr .

For sufficiently fast rotators, the thermal skin depth can be exceedingly small. For the special case of δ ! D, the macroscopic stresses are anticipated to asymptotically ap-Fig. 10: The non-monotonic dependence of maximum macroscopic surface stress on the rotation period is shown for three different rock sizes on an S-type asteroid at 1 AU. For each rock size, the particular rotation period that results in the highest macroscopic surface stress is termed the critical rotation period, Pcr, which is associated with a corresponding critical thermal skin depth δcr. From these calculations, it is apparent that δcr " D{2.

proach Σ8

s " βαE∆ Ts {16p1 ´νq as implied by Eq. ( 20). In this regime, the macroscopic thermal stresses are anticipated to be size independent. Interestingly, the numerical calculations for Σs shown in Fig. 10 for various rock sizes seem to be converging, i.e. becoming size independent, for extremely fast rotators (P À 10 ´1 ´10 0 hours). The larger the rock sizes, the easier it is to realistically achieve the condition δ ! D. As such, we anticipate the peak thermal stresses in rocks that are 10s of meters in size and larger will tend to be nearly size-independent for planetary bodies with rotation periods faster than about an Earth day. On the Moon, the thermal stresses are anticipated to become sizeindependent for rocks larger than about 100 meters. Lastly, since thermal stresses in this size-independent regime scales only with thermoelastic material properties and the diurnal surface temperature variation, i.e. Σ8 s 9∆ Ts , the scaling of diurnal surface temperature variation with spin rate will fully govern the spin rate dependence of Σs , which itself was shown in Fig. 2 to scale roughly with ? P for sufficiently fast rotators and independent of rotation period for sufficiently slow rotators (such that ∆ Ts Ñ T eq q.

The scaling relations derived in section 5.3 have proven invaluable in providing fundamental insights into the nonmonontonic dependence of macroscopic thermal stresses first reported in [START_REF] Ravaji | The effect of rotation period on thermal stress weathering[END_REF] and interpretation of Fig. 10. For example, it is now clear that there is a critical spin rate for each rock size for which macroscopic thermal surface stresses will be maximum. For spin rates faster than this critical spin rate, the macroscopic thermal stresses are governed primarily by thermoelastic properties, ∆ T , and δ{D with this size-dependence vanishing for sufficiently large rocks and/or sufficiently fast rotators. For spin rates slower than the critical spin rate, the macroscopic thermal stresses are governed primarily by thermoelastic properties, |∇ 2

Ts |, and the square of the rock size, i.e. D 2 . Beyond these fundamental insights, the scaling relations are useful approximations that are remarkably accurate when compared to full numerical calculations as well as three dimensional finite element calculations. For example, the dotted lines in Fig. 10 are predictions of the scaling relation expressed in Eq. ( 20), which exhibits remarkable agreement with the numerical calculations (markers in Fig. 10). For completeness and comparison, the dashed and solid lines in Fig. 10 make use of higher-fidelity analytic approximations, i.e. first-order and second-order Fourier series approximations derived in Appendix A. The first-order and second-order Fourier series exhibit better agreement with the numerical calculations, but are a bit more cumbersome than Eq. ( 20). The sharp transitions in the second-order Fourier series are associated with a transition in the time of day when peak macroscopic thermal surface stresses are achieved.

The combined effect of spin rate and orbital distance on the maximum microscopic thermal stresses are shown as a contour map in Fig. 11 for a 10 cm boulder on an S-type asteroid. The contour map is generated via combination of Eq. ( 20) and Eq. ( 11). As anticipated Σs is greatest for asteroids that pass close to the sun and decays monotonically with increasing orbital distance in a manner similar to that of ∆ Ts , i.e. as r ´2 for fast rotators and r ´1{2 for slow rotators according to Eqs. ( 9)-( 11). [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF] reported similar scaling dependencies of r ´3{2 and r ´1.2 on Σs for rotation periods of 1 and 5 Earth days, respectively, and slightly stronger dependencies for faster rotating bodies. [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF] provided a map of finite element calculations of Σs similar to Fig. 11, but did not observe the non-monotonic dependence of Σs on spin rate. The reason for that is that [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF] calculations are for a much larger boulder size, D " 1 m (compared to D " 10 cm in Fig. 11). Our results and the calculations of [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF] agree for equivalent parameters. As discussed earlier, the critical spin rate where the spin rate dependence of Σs reverses corresponds to spin rate required to produce a thermal wave that penetrates halfway through the rock. For the 10 cm rock modeled in Fig. 11, this condition is met for a rotation period of about 0.1 Earth days. However, for the 1 m rock modeled in [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF] a rotation period of about 12 Earth days is required to meet the condition, and the map shown in [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF] does not extend to such slow spin rates. Lastly, it is interesting to note that there are no known asteroids in the regions of Fig. 11 with Σs Á 4 MPa, i.e. r À 0.75 AU and 10 ´2 Earth days À P À 0.5 Earth days. It is possible that thermal stress weathering on centimeter-and meter-sized is highly effective on bodies with such rotation periods and orbital distances.

Summary and conclusions

The efficacy of thermal stress weathering is driven by both macroscopic (rock-scale) and microscopic (grain-scale) thermal stresses. Here we carried out extensive numerical calculations and analytic analysis of the [START_REF] Delbo | Thermal fatigue as the origin of regolith on small asteroids[END_REF] thermal stress weathering model. Microscopic thermal stresses were found to be governed primarily by thermoelastic material properties, the volume fraction of heterogeneities, the contrast between thermal expansion coefficients of heterogeneities, and the maximum diurnal surface temperature variation as was also stated by [START_REF] Molaro | Grain-scale thermoelastic stresses and spatiotemporal temperature gradients on airless bodies, implications for rock breakdown[END_REF] and [START_REF] Eppes | Mechanical weathering and rock erosion by climate-dependent subcritical cracking[END_REF]. Macroscopic thermal stresses are more complex and were found to be governed primarily by thermoelastic material properties, the ratio of thermal skin depth to rock diameter, the maximum diurnal surface temperature variation, and the second-order spatial derivatives of the temperature field, with many of these factors also being articulated in [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF]. Beyond these fundamental insights, here we have derived a number of scaling relations that are practical approximations of the magnitude of both microscopic and macroscopic thermal stresses. These scaling relations are remarkably accurate when compared to both the numerical calculations as well as full finite element calculations [START_REF] Molaro | Grain-scale thermoelastic stresses and spatiotemporal temperature gradients on airless bodies, implications for rock breakdown[END_REF][START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF]. In addition, these scaling relations accurately capture the effect of both thermoelastic and thermophysical properties on microscopic and macroscopic thermal stresses.

The peak microscopic (grain-scale) thermal stresses are reported to be roughly 2-4 MPa on many near-Earth asteroids and about 36 MPa on the Moon, which is in good agreement with two-dimensional mesoscale finite element calculations of [START_REF] Molaro | Grain-scale thermoelastic stresses and spatiotemporal temperature gradients on airless bodies, implications for rock breakdown[END_REF] (after necessary modifications are invoked to account for the different boundary conditions utilized). Likewise, the magnitude of peak macroscopic thermal stresses on the Moon is reported here to respectively be about 5 MPa and 8 MPa on a D " 0.5 m and a D " 0.7 m rock, which is in remarkably good agreement with three-dimensional finite element calculation of [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF]. Stress magnitudes reported in this paper through Fig. 6 to Fig. 11 are generally subcritical compared to the tensile strength of rocks [START_REF] Eppes | Mechanical weathering and rock erosion by climate-dependent subcritical cracking[END_REF][START_REF] Ravaji | Introducing optimized validated meshing system for wellbore stability analysis using 3D finite element method[END_REF], but large enough to trigger fatigue process. That said, thermal fatigue lifetime calculations were beyond the scope of the current work.

One important contribution of this paper has been the classification of a body's spin rate based on the thermal parameter, Θ, where different fundamental physics and scalings dominate, i.e.

• Thermally fast rotators exhibit Θ " 1

• Thermally slow rotators exhibit Θ ! 1

• Intermediate rotators exhibit Θ " 1.

For sufficiently fast rotators with Θ Á 10, e.g. P À 6 hr at 1 AU from the Sun, the diurnal surface temperature variation, and thus, the peak microscopic thermal stresses are inversely proportional to the thermal parameter Θ, and roughly constant for constant values of Θ. As such, microscopic thermal stresses will increase roughly with the square root of the rotation period ? P and decay with orbital distance as r ´2. For sufficiently slow rotators with Θ ! 1, the peak microscopic thermal stresses are roughly independent of spin rate and are expected to decay with orbital distance as r ´1{2 . In addition, the two regimes exhibit different dependencies on thermophysical parameters, e.g. the diurnal surface temperature variation and peak microscopic thermal stresses are inversely proportional to thermal inertia Γ and independent of emissivity on thermally fast rotators, but are independent of thermal inertia Γ and dependent on emmisivity on thermally slow rotators. These dependencies on thermophysical properties and orbital parameters can be seen clearly via combination of the scaling laws shown in Eq. ( 11) and Eq. ( 13), which provides an attempt to generalize observations of [START_REF] Molaro | Grain-scale thermoelastic stresses and spatiotemporal temperature gradients on airless bodies, implications for rock breakdown[END_REF], i.e.

• On thermally slow rotators with Θ ! 1: ∆ Ts " 2 3 T eq " 2 3

4 c 1 ´A σ B S @ r 2 σT M s " 2 3 K ˚f pv i q∆α 4 c 1 ´A σ B S @ r 2
• On thermally fast rotators with Θ " 1:

∆ Ts " 1.12 ˆTeq Θ " 1.12 p1 ´AqS @ Γ ? ωr 2 σT M s " 1.12 ˆK˚f pv i q∆α p1 ´AqS @ Γ ? ωr 2 .

An additional classification introduced here regards the relationship between a structural length scale, the rock diameter D, and a thermal length scale, the thermal skin depth δ. Here, large rocks in a thermal stress context are those whose structural length scale is larger than the thermal length scale, i.e. D ą δ, and conversely small rocks are those with D ă δ. For sufficiently small rocks, i.e. D ă δ, the peak macroscopic (rock-scale) thermal stresses are approximately proportional to the magnitude of the second-order gradient in surface temperature, i.e. ∇ 2

Ts , as shown in Eq. ( 19), which scales roughly with the inverse square root of the rotation period (1{ ? P ) on fast rotators (Θ ą 1) and as the inverse of rotation period on (1{P ) on slow rotators (Θ ă 1), as discussed in subsection 4.3. This difference in the dependence on ∇ 2

Ts instead of ∆ Ts is the most crucial difference between macroscopic thermal stresses on small rocks and microscopic thermal stresses, and leads to opposite dependencies on spin rate. The macroscopic thermal stresses on small rocks are strongly dependent on rock size, and scale with the square of the rock size as shown in their scaling relations resulting from Eq. ( 19) combined with Eq. (11-12):

• Small rocks, i.e. D À δ, on slow rotators, i.e. Θ ! 1:

Σs " 1 27

αE 1 ´ν ˆD δ ˙2 4 c 1 ´A σ B S @ r 2
• Small rocks, i.e. D À δ, on fast rotators, i.e. Θ " 1:

Σs " 1.12 ˆ1 18 αE 1 ´ν ˆD δ ˙2 p1 ´AqS @ Γ ? ωr 2 .

On the other hand, for sufficiently large rocks, i.e. D " δ, the peak macroscopic (rock-scale) thermal stresses are observed to follow fairly similar trends as those discussed for peak microscopic (grain-scale) thermal stresses, since both are primarily drive by the diurnal surface temperature variation ∆ Ts , with the scaling laws for both exhibiting similar functional forms on both fast and slow rotators as seen by combination of the scaling laws shown in Eq. ( 11) and Eq. ( 21), which provides an attempt to generalize observations of [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF], i.e.

• Large rocks, i.e. D Á δ, on slow rotators, i.e. Θ ! 1:

Σs " αE 1 ´ν f ˆδ D ˙4 c 1 ´A σ B S @ r 2
• Large rocks, i.e. D Á δ, on fast rotators, i.e. Θ " 1:

Σs " 1.12 ˆ3 2 αE 1 ´ν f ˆδ D ˙p1 ´AqS @ Γ ? ωr 2 .

The primary difference between the peak macroscopic (rockscale) and microscopic (grain-scale) thermal stresses in large rocks (D Á δ) is that the macroscopic thermal stresses depend on δ{D, which itself depends on both the spin rate ω and the thermal diffusivity. For extremely large rocks D " δ the macroscopic thermal stress becomes insensitive to δ{D, i.e. f `δ D ˘Ñ 1, and the dependence of peak macroscopic and microscopic thermal stresses on spin rate and thermophysical properties become essentially identical. It is also worth noting that the dependence of both microscopic thermal stresses and macroscopic thermal stresses on orbital distance is essentially identical. For all cases, we report that the thermal stresses should decay as r ´2 on very fast rotators and as r ´1{2 on very slow rotators. Intermediate rotators will have scaling dependencies that fall between these two limiting cases. This is consistent with calculations of [START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF] who reported an orbital distance scaling dependence with an exponent that weakened with increasing rotation period, i.e. r ´1.65 , r ´1.5 , and r ´1.2 for a 1 m boulder on a body with a period of 0.1, 1, and 5 Earth days, respectively.

An additional important contribution of this paper is the physical rationalization of a non-monotonic dependence of macroscopic thermal stresses on spin rate, first reported by [START_REF] Ravaji | The effect of rotation period on thermal stress weathering[END_REF]. Here we postulate and confirm that there is a critical spin rate where macroscopic thermal stresses will be greatest for a given rock size. This critical spin rate corresponds to the spin rate necessary to induce a thermal wave that penetrates halfway through the rock, i.e. δ cr " D{2. Spin rates either faster or slower than the critical rate will result in smaller macroscopic thermal stresses. In the above derived scaling law, this non-monotonic dependence arises through the competition of 1{ ? ω and f `δ D ˘. As the spin rate slows, the thermal skin depth δ increases, which results in a monotonic decrease of f `δ D ˘, whereas 1{

? ω monotonically increases as the spin rate slows. The product of these two monotonic, but opposite trending functions captures the non-monotonic dependence on spin rate.

Acknowledgments

MD acknowledges support from the Centre National d Études Spatiales, as well as the Academies of Excellence on Complex Systems and Space, Environment, Risk and Resilience of the Initiative d'EXcellence "Joint, Excellent, and Dynamic Initiative" (IDEX JEDI) of the Université Côte d'Azur. We note that there are no data sharing issues since all of the numerical information is provided in the figures produced by solving the equations in the paper.

Appendix A. Fourier series analysis

One convenient, general harmonic solution of Eq. ( 7) for the one-dimensional spatiotemporal temperature field is T pz, tq " Ts `1 2 ∆ Ts

Fig. 3 :

 3 Fig.3: Diurnal surface temperature variation ∆ Ts (normalized by equilibrium model temperature Teq) as a function of the thermal parameter Θ. There is remarkable agreement between scaling law provided by Eq. (11) shown as solid red line, and numerical data shown as blue circles. The numerical data includes predictions at various orbital distances and for 2 sets of thermophysical material properties, i.e. ordinary and carbonaceous chondrites.

Fig. 4 :

 4 Fig. 4: Normalized ∆ Ts, |∇ Ts| and |∇ 2 Ts| at the surface of an Stype asteroid versus orbital distance and thermal parameter Θ. The numerical calculations are normalized by the corresponding value at 1 AU.Here the rotation period is taken to be 6 hours.

Fig. 6 :

 6 Fig. 6: Contour plot demonstrating the combined effect of rotation rate and orbital distance (assuming zero eccentricity and inclination) on the maximum microscopic thermal stress, i.e. σT M s , within a 10 cm rock on an S-type asteroid. For reference, the rotation rate and average orbital distances of various planetary bodies (Hergenrother and Whiteley, 2011) is shown.

  (b) with Fig. 7 (f).

Fig. 7 :

 7 Fig. 7: Time evolution of (a)-(b) macroscopic thermal stress, (c)-(d) temperature, and (e)-(f) time-rate-of-change of temperature at various depths for a 10 cm boulder on an S-type asteroid rotating with a period of P " 0.6 hours (Θ " 30) and P " 103 hours (Θ " 2.4).

Fig. 8 :

 8 Fig.8: Comparing 3D finite element calculations of[START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF] for a boulder on the Moon surface with analytical approaches provided here via Eq. A.3 and Eq. (21) for the maximum stress at the surface, Σs, and Eq. (22) for half rock diameter, ΣD{2 , with material properties adopted from[START_REF] Molaro | Thermally induced stresses in boulders on airless body surfaces, and implications for rock breakdown[END_REF]. Although stress magnitudes are different (most likely because of the difference between 1D and 3D conditions), trend of the plots are in fairly good agreement.

Fig. 11 :

 11 Fig. 11: Contour plot demonstrating the combined effect of rotation rate and orbital distance (assuming zero eccentricity and inclination) on the maximum macroscopic thermal stress, i.e. Σs, within a 10 cm rock on an S-type asteroid. For reference, the rotation rate and average orbital distances of various planetary bodies (Hergenrother and Whiteley, 2011) is shown.

  P ´1 ş P {2 ´P {2 T s ptqdt denotes the mean surface temperature through a diurnal cycle. In general, the Fourier series coefficients utilized in Eq. A.1 are computed as A.1 in Eq. (3) leads to an analytic equation for the temporal variation of the macroscopic thermal stress on the surface, i.e.

  Ts {dt| is easier to measure experimentally than |∇ 2

		Ts |;
	however, |∇ 2	Ts | is a more fundamental factor governing
	macroscopic stresses. Nevertheless, in subsequent sections,
	|∇ 2