Uncertainties in physical properties of Itokawa-like asteroids widen constraints on their formation time

One of the outstanding questions in planetary science is to determine how the fundamental mechanical and physical properties of materials determine the thermal evolution of asteroids, and which properties have the greatest influence. We investigate the effects of uncertainty in the material properties of asteroid parent bodies on the ability of thermal evolution models to constrain the sizes and formation times of ordinary chondrite parent asteroids. A simple model is formulated for the thermal evolution of the parent body of asteroid 25143 Itokawa. The effects of the uncertainties in the values specified for specific heat capacity, thermal diffusivity, and aluminum abundance are determined. The uncertainties in specific heat capacity and aluminum abundance, or heat production more generally, are found to both have significant and approximately equal effects on these results, substantially widening the range of possible formation times of Itokawa’s parent body. We show that Itokawa’s parent body could have formed between 1.6 and 2.5 million years after the origin of calcium–aluminum inclusions with a radius larger than 19 km, and it could have formed as early as 1.4 millions years, as late as 3.5 million years, or with a radius as small at 17 km if more lenient definitions of uncertainty in aluminum abundance are considered. These results stress the importance of precise data required of the material properties of a suite of LL type 4-6 ordinary chondrite meteorites to place better constraints on the thermal history of Itokawa’s parent body.


Introduction
The current near-earth asteroid and main-belt smallbody populations include both leftover rocky and ice-rich bodies formed in the early solar system and processed, differentiated fragments of large asteroids. A majority of the chondritic materials were formed and partially processed approximately < 10 million years (Ma) after time 0 (defined as the time when calcium-aluminum inclusions, CAIs, originated) (Desch et al. (2018) and references therein). The physical and chemical characteristics of the chondritic materials that make up the asteroids are a function of precursor material compositions and subsequent alteration mechanisms that occurred after the asteroids grew to a particular size. The primary alteration mechanism for >10-km-sized bodies is due to the presence and decay of 26 Al, although there is conclusive evidence that impact and shock processing were prevalent in the early solar system [McSween et al. (2002) and references therein]. Thermal metamorphism of ordinary chondrites (OCs) occurred during the first 5 to 50 Ma after a parent body's accretion, depending on the size of the body (Huss et al. 2006).
Theoretical and computational models that calculate the thermal evolution experienced by an early-formed solar system body require experimental measurements of asteroid material to constrain evolutionary details about the asteroid including the body's size and time of formation. A simulated thermal history of a parent body includes several free or unknown parameters, and conducting many simulations can determine what combinations of free parameters result in thermal evolutions that comply within the material constraints and observations. For OCs some of these studies include Henke et al. (2012a), Wakita et al. (2014), Sugiura and Fujiya (2014), Doyle et al. (2015), Blackburn et al. (2017), and Gail and Trieloff (2019).
The ability of this approach to make significant predictions about asteroid parent bodies is contingent on the accuracy and precision of thermal constraints determined from meteorites in our inventory and the completeness of the thermal evolution model and its fixed parameters. Even for the simplest models and types of parent bodies, there are many possible sources of significant error. These include omission or oversimplification of significant physical phenomena that impact the body's evolution, errors arising from the approximate solutions of model equations and other computational techniques, and uncertainties in the values and equations used to describe the fixed parameters of the model. The studies that utilize these models typically examine some of these sources of error, but the combined quantitative effect of these errors on their results can be prohibitively difficult to determine. Previous discussions of error in parent body thermal evolution models include bounding the error of the finite difference methods by comparison to analytically solved bodies such as in Wakita et al. (2014), an approximate accounting of the effects of uncertainties in peak central temperatures and aluminum abundance on formation time ranges in Sugiura and Fujiya (2014), and various qualitative remarks on the magnitude of certain errors in other studies. This study focuses on characterizing the effects of uncertainty in material properties on the thermal histories and predictions made by such models.
We will begin by motivating and reconstructing the thermal evolution model of the parent body of asteroid 25143 Itokawa as described in Wakita et al. (2014). The constraints on the formation time and size of Itokawa's parent body predicted by this model will then be recalculated while accounting for the uncertainties in the values measured for specific heat capacity, thermal diffusivity, and aluminum abundance in OCs. The model of Itokawa discussed here is relevant to many of the recent and ongoing studies of Itokawa's characteristics and origins, such as the nature of organics and water in the asteroid (Jin and Bose 2019;Chan et al. 2021). The results of this study are also applicable to other undifferentiated asteroids that are of current interest to the community, e.g., C-type asteroids Bennu and Ryugu.

Thermal evolution model and parameters
We based our model and approach on Wakita et al. (2014), which investigated the parent body of asteroid 25143 Itokawa. The majority of the details and model parameters are assumed from their paper. The only notable change is our use of an implicit finite difference approximation to solve the heat conduction equation. The implicit finite difference approximation was chosen because of the improved stability and convergence of its solutions, and has been used and discussed in other thermal evolution models (Henke et al. 2012a).
The body is assumed to be spherically symmetric, made entirely of a homogeneous material, and only heated by the decay of radionuclide 26 Al with a short half-life (halflife = 720,000 years; Miyamoto et al. (1981)). The use of additional radionuclide heat sources, such as 60 Fe or 40 K (Henke et al. 2012a;Gail and Trieloff 2019), were also briefly considered, but testing their inclusion in a few simulations showed they were not significant for this model and study, as also found in Wakita et al. (2014).
The material of the body is also assumed to be nonporous and incompressible such that the body has a constant size over time. A significant degree of porosity in the material of asteroid parent bodies can have notable effects on their evolution, as is discussed in other thermal evolution models. Relevant examples include the model for the parent body of C-type asteroid Ryugu by Neumann et al. (2021), which accounts for initial parent body porosities of 50% and higher, and the L chondrite parent body model by Gail and Trieloff (2019), which uses a value of 25% for porosity. Parent body porosity is primarily determined by the microporosity of relevant samples, and particles from Itokawa have been noted to have an unusually small microporosity, an average of only 1.4% from the analysis of Tsuchiyama et al. (2011). For this reason, the assumption of zero porosity for Itokawa's parent body is unlikely to significantly impact the results of this thermal evolution model.
The specific heat capacity of the material is given by the equation (Yomogida and Matsui 1984): where c p is specific heat capacity and T is the temperature. The thermal diffusivity or thermometric conductivity, κ , which has the relation to other material properties κ = K ρc p where K is thermal conductivity and ρ is density, is given for this material by the equation (Yomogida and Matsui 1984): The spherically symmetric heat conduction equation, which describes how temperature changes over space and time in the body, is Using thermal diffusivity (Eq. 2) instead of thermal conductivity and assuming that the density stays constant allows us to cancel out the dependence on density from each term, making the results of this model independent of the value used for material density and simplifying the equation to: To solve the heat conduction equation, we followed the approach of Henke et al. (2012a) and used a fully implicit Crank-Nicolson method to approximate thermal conduction over time. We imposed a Neumann boundary condition at the center of the body, which disallows any heat transfer through the center. The surface of the body had a Dirichlet boundary condition, being held at a constant temperature of 200 K as in Wakita et al. (2014). The approximation solved for 300 nodes spaced evenly across the body's radius, as using more nodes or unequal node spacing did not have a significant effect on the results. The timestep varied between a minimum of 10 years and a maximum of 1000 years. We shortened the timestep to its minimum during periods of fast thermal changes, when the temperature of any node changed by more than 0.3% in a single timestep, and increased the timestep to its maximum during periods of slow thermal changes, when the temperature of every node changed by less than 0.01% in a single timestep (Henke et al. 2012a).
where the parameter A is the magnitude of the heat energy of the decay and the base value used is 2.034 × 10 −7 W/kg, derived from the values used in the Wakita et al. (2014) model. Specifically, the heat production can be calculated from A = X m f E τ (Henke et al. 2012a) where X is the mass fraction of aluminum at 1.18 wt.%, f is the 26 Al/ 27 Al ratio of 5 × 10 −5 at the time of CAI formation, E is the thermal energy of 3.16 MeV per decay of an atom of 26 Al, τ is the time constant of the decay at 1.039 Ma from the half-life of 26 Al of 0.72 Ma, and m is the atomic mass of aluminum at 26.98 u.

Reported uncertainties in physical properties and our assumptions
While there are many possible sources of error in the model described above, this study focuses on understanding the effects of the uncertainties in the material parameters of the model, primarily the values and equations used for thermal diffusivity, specific heat capacity, and aluminum abundance.
First, reasonable uncertainty ranges for the base values and equations of these material properties must be established. Because the base values for specific heat capacity and thermal diffusivity are functions of temperature, Eqs.
(1) and (2), respectively, uncertainty in these properties could be defined in several ways. For the purposes of this study, the uncertainty in these properties will be implemented by multiplying each term in that property's temperature-dependent expression by the same factor, which is determined from the uncertainty range. Other implementations of uncertainty, which could include scaling the different terms of the expressions by unequal factors or altering the order or form of the temperature dependence, are not investigated. For example, a ± 10% uncertainty in a property and equation would be implemented and discussed as a range of multiplicative factors from 0.9 to 1.1. The uncertainty in aluminum abundance will also be referred to as a range of multiplicative factors, in reference to its base value, for consistency.
Ideally, the uncertainty of a model parameter would be explored with simulations that vary the parameter by many different points in its range of multiplicative factors. However, our simulations showed that it was only the extreme, end-member values of these ranges that produced the most extreme changes in the results, and the curves and numbers which characterize these results changed smoothly and monotonically as the uncertain parameters were varied, so the majority of our simulations only used either base or end-member values for these uncertain parameters. Overall, about 1900 simulations were performed. For combinations of parameters at their base value or end-member uncertainty value (data presented explicitly here), 80-120 simulations were performed per unique combination of uncertain parameters in order to sufficiently span the free parameter space of radius and formation time values. For combinations of uncertain parameters that used only intermediate values of their uncertainty ranges, only ∼ 20 simulations were performed to confirm that the behavior of these combinations agree with the overall trends.

Specific heat capacity
Temperature-dependent specific heat capacity equations are important to use in thermal evolution models because of the large degree in which the value changes over the temperatures relevant to parent body modeling (Ghosh and McSween Jr. 1999). Several papers use a constant value for specific heat capacity for some or all of their simulated body's material [(e.g., Blackburn et al. (2017), Doyle et al. (2015), Fujiya et al. (2012)]. Those that use a temperature-dependent value typically use the experimentally fitted equation from Yomogida and Matsui (1984) and Yomogida and Matsui (1983), as in Wakita et al. (2014) and our model, or use weighted sum methods to approximate the heat capacity at a given temperature from the heat capacities of its base components at that temperature, as in Henke et al. (2012a).
With the equation from Yomogida and Matsui (1984) (Eq. 1) as our base value, we chose a multiplicative factor range from 0.9 to 1.1 for specific heat capacity (Table 1). The primary justification for this range is a comparison of the base equation to other models and measurements of OC specific heat capacity at various temperatures. The specific heat capacities for H and L chondrites used in Henke et al. (2012a) go both above and below the base equation by approximately 10% at certain temperatures between 1000 and 1300 K, and measurements of the specific heat capacity of LL and L chondrites shown in Flynn et al. (2018) also go above and below the base equation by ∼ 10% at temperatures between 300 and 900 K. Additionally, the data that are the basis of the Yomogida and Matsui (1984) equation are primarily not from LL chondrite samples and do not include measurements from material temperatures above 500 K, both of which justify a reasonable amount of uncertainty in the use of the equation in this particular model.

Thermal diffusivity
We used the thermal diffusivity equation from Yomogida and Matsui (1983), (Eq. 2), specifically the fit to the LL6 OC sample of St. Severin as used in Wakita et al. (2014). The thermal diffusivity measurements reported in Yomogida and Matsui (1983) are stated to have experimental uncertainties of 10%, but there are additional uncertainties which we took into account. First, their measurements only go up to 500 K, and undifferentiated OC bodies can reach temperatures up to 1300 K. Second, their other LL6 sample, St. Lawrence, has fit parameters that differ from St. Severin by more than 10%. Lastly, additional studies have shown that the thermal conductivity of pristine material in asteroid parent bodies could be multiple times greater than the measured conductivity of analogous meteorite samples, due to the considerable extent to which shock modification of meteorites that occurred after the material was part of its parent body reduces the thermal conductivity of that material (Gail and Trieloff , 2018;Henke et al. , 2016).
An independent value for thermal conductivity is not explicitly given in this study, but rather it can be determined by the relation K = κρc p . Thermal diffusivity and specific heat capacity are given by their temperaturedependent expressions, Eqs. (1) and (2), and the value of density used is 3400 kg/m 3 (Tsuchiyama et al. 2011), as in Wakita et al. (2014). Using these base values, a conductivity of 1.53 W/mK is determined for the material at room temperature (300 K), and a value of 1.68 W/mK if the largest value of specific heat capacity as defined in "Specific heat capacity" is used. The thermal conductivity model of Henke et al. (2016) finds that the room temperature conductivity of pristine ordinary chondrite material with vanishing porosity would be 4.89 W/mK for H chondrites, 4.20 W/mK for L chondrites, and 3.78 W/mK for LL chondrites, which we take as upper bounds for the derived thermal conductivity of our material. Because the uncertainty in material density is not considered here, and the uncertainty range for specific heat capacity defined in "Specific heat capacity" is already at the extremes of what is reasonable for ordinary chondrite material, this upper bound of the thermal conductivity value can therefore be related directly to the upper bound of the value for thermal diffusivity. Comparing the value of 1.68 W/mK for our material to these upper bounds of thermal conductivity, which are up to 2.9 times greater, and taking into account the additionally mentioned sources of error for thermal diffusivity justifies a multiplicative variation range of 0.9 to 3.0 for thermal diffusivity.

Aluminum abundance
The value of the heat generation parameter A in the heat production equation (Eq. 5), which is the magnitude of the heat power released by the decay of 26 Al, is directly proportional to the abundance of aluminum in the asteroid material. Wakita et al. (2014) used a value of 1.18 wt.% based on the measurements of Jarosewich (1990), which measured the aluminum abundance in LL6 samples to range from 1.04 to 3.73 wt.%. We used the same 1.18 wt.% as our base value, and comparing to that range of measurements gives a multiplicative variation range of 0.88-3.16 (Table 1). However, this measurement of 3.73 wt.% aluminum comes from a single measurement of the St. Mesmin meteorite and is an outlier compared to the other OCs measured in that dataset (Jarosewich 1990). Because of this, a more conservative uncertainty range was also considered. Sugiura and Fujiya (2014) states that aluminum abundance of chondrites can vary up to 15% from their value of 1.19 wt.%, which is generally consistent with the majority of the relevant Jarosewich (1990) measurements, and so for our base value this would give a smaller range of variation from 0.86 to 1.16 (Table 1). Heat capacity, thermal diffusivity, Al abundance (Jarosewich 1990) 0.9-1.1, 0.9-3.0, 0.88-3.16 a 1.64-3.54 16.5 Heat capacity, thermal diffusivity, Al abundance  0.9-1.1, 0.9-3.0, 0.71-0.87 a 1.43-2.20 19.6 Nakamura et al. (2014) analyzed 48 Itokawa particles and reported that the bulk composition of the particles was 0.937 ± 0.095 wt.% aluminum. This is less than the aluminum abundances for all of the LL chondrite samples reported in Jarosewich (1990), and less than all but two L chondrite samples reported therein. Because this finding is a significant outlier compared to analogous materials and was determined from a much smaller amount of material (< 1 mm) than in most of the measurements reported in Jarosewich (1990), we report our model results for this range, 0.71-0.87 times our base value (Table 1), separately.
It should also be noted that a change in the aluminum abundance of the body's material would also have some effect on the specific heat capacity and thermal diffusivity of the material, as would any change in the material's composition. However, because the aluminum abundance is only a small portion of the material's composition the change in specific heat capacity and thermal diffusivity due to varying the aluminum abundance is likely small and can be safely ignored.

Results
Following the criteria of Wakita et al. (2014), we considered a simulated body to be a suitable parent body for asteroid Itokawa if the center of the body, which is the highest temperature point at all times, fulfilled three requirements or constraints: its peak temperature was greater than 1073 K, its peak temperature was less than 1273 K, and its temperature at 7.6 Ma after the formation of CAIs was greater than 973 K.
The only free parameters left in this thermal evolution model are the radius and formation time of the simulated body, and only certain combinations of values for these two parameters result in suitable parent bodies of Itokawa. These combinations of values can be analyzed by finding the points in the two-dimensional parameter space which result in simulated bodies that fulfill all of the above criteria. As in Wakita et al. (2014), this is accomplished by first finding the peak central temperatures over this parameter space and constructing the isothermal lines for 1073 and 1273 K, then finding the central temperatures at 7.6 Ma over the parameter space and constructing the isothermal line for 973 K, and finally overlaying the three isotherms. The area enclosed by the three lines is then all of the possible combinations of radius and formation time values that fulfill the criteria. Note that for this and many thermal evolution models, the peak central temperatures of simulated bodies will cease to be affected by a change in body radius for radius values greater than a certain value. For most simulated bodies considered here, the peak central temperature of a 50-km radius body is usually within 1% of the peak temperature of a larger but otherwise similar body. For the bodies we consider with large thermal diffusivity values, however, it is not until about 70 km that increasing the radius ceases to significantly affect the peak central temperature. Because the 1073 and 1273 K isotherms of peak central temperature are non-intersecting and parallel above radius values of about 70 km, the three thermal constraints together cannot bound a finite area, but instead the 973 K isotherm of temperature at 7.6 Ma serves to mark the lower boundary of the area and there is no constraint to mark the upper boundary. For this reason, the three thermal constraints are unable to bound a maximum radius value for the parent body of Itokawa, but are able to determine a minimum radius and earliest formation time by the intersection of the 1273 K peak temperature isothermal line and the 973 K at 7.6 Ma isothermal line, and a latest formation time by the 1073 K peak temperature isothermal line.
For every combination of end-member parameter variations, 80-120 simulations were used to span the free parameter space and construct the three isotherms which give the formation time range and minimum radius for that parameter combination. Accuracy of these isotherms in predicting the model behavior is critical to the discussion on the effects of the parameter variation, and later additional simulations have confirmed that the formation time ranges and minimum radii stated throughout and in Table 1 are true to the model on the order of 0.01 Ma and 0.1 km. The accuracy of these isotherms was improved by doing two sets of simulations for each parameter combination: 40-60 simulations spanned the free parameter space of interest, and then preliminary isotherms were constructed from that data, and finally 40-60 additional simulations where performed with free parameter values much closer to the boundaries of those preliminary isotherms. We report only two significant figures when stating our findings for the parent body of Itokawa.

Testing the model for Itokawa
The initial simulations (120 in total) used only the base values and equations for all the model parameters and simulated bodies with a variety of radius and formation time values shown by the positions of the gray dots in Fig. 1. The peak central temperatures and central temperatures at 7.6 Ma were interpolated over the two-dimensional parameter space of interest from these points, and then the three isothermal lines were drawn to locate their intersections and bounded area.
While there are visible interpolation defects in the resulting isotherms, which are not as smooth as would be expected, the magnitude of these defects is small and therefore not considered to be a major source of error in the results. The area between all three isotherms shows the range of radius and formation time values possible for Itokawa's parent body as predicted by our base model without incorporating any uncertainty in material properties, and is given as our result with no parameter variation in Table 1 Figure 10 of Wakita et al. (2014). This was an important first step to ensure that our Itokawa model agreed with published work, and confirmed the general accuracy of our finite difference solver and our use of interpolation.
Our finite difference solver was also incorporated into a model that follows the parameters and design of Miyamoto et al. (1981), whose simplified parent bodies have analytically solved thermal evolutions. Our finite difference solution was able to reproduce the temperatures throughout the body to within 1% of the analytical solution for all except the earliest stages of the evolution. The evolution immediately following formation is not as accurate due to our choice in timestep. This does not significantly affect the accuracy of the peak temperatures and temperatures after several million years, however, which are the primary concern for this study.

Simulations and isotherms including uncertainties
The effects of varying only specific heat capacity within its uncertainty range were found by following the same steps as above except with modified specific heat capacity equations. Two additional sets of 90 simulations were conducted with specific heat capacity multiplied by the two end-member values of its multiplicative factor range. The results are shown in Fig. 2 and stated in Table 1. As seen in Fig. 2a, varying specific heat capacity has only a minor effect on the minimum radius size of Itokawa's parent body, with the combined isothermal lines of the varied simulations shifting and stretching only slightly vertically. There is a much more significant effect on the range of valid formation times, seen by the large shifts towards earlier or later times by the bounded area of the isothermal lines corresponding to higher or lower heat capacity, respectively. With only specific heat capacity varied, the possible formation times can be as early as 1.77 Ma and as late as 2.35 Ma. Though the relationship between the amount of uncertainty and the change in results is not necessarily linear, the ± 10% uncertainty in specific heat capacity resulted in an additional ∼ ± 0.1 Ma added onto our formation time range, which is a decrease in the lower bound and an increase in the upper bound by ∼ 5% each. Additionally, this amount of uncertainty in specific heat capacity lowered the minimum radius by approximately 2.5%.
After the creation of these combined isothermal lines, 2 additional simulations with end-member specific heat capacity values and formation times on the boundaries of this predicted range were conducted. The central temperatures of these simulated bodies are shown as a  Wakita et al. (2014) with no uncertainties in model parameters taken into account. The gray dots are locations in the free parameter space at which simulations were performed. The peak central temperatures achieved in these simulations were linearly interpolated to create the color map pictured. Two of the isotherms shown correspond to 1073 and 1273 K on the peak central temperature map shown in color, and the third isotherm corresponds to 973 K on a map of central temperatures at 7.6 Ma which is not pictured function of time in Fig. 2b, c, and the radius and formation time values are marked by diamonds in Fig. 2a. Figure 2b shows the thermal evolution of a 21-km radius body formed at 1.78 Ma with a specific heat capacity at 1.1 times the base equation, and Fig. 2c shows a 50-km radius body formed at 2.34 Ma with a specific heat capacity at 0.9 times the base equation. Both of these simulations fulfill all three criteria, but are also clearly close to failing at least one of the criteria. This is expected for bodies with end-member formation time values and demonstrates that the isothermal lines drawn are consistent with our simulation results and not significantly impacted by the interpolations used in creating them.
Varying thermal diffusivity was performed in the same manner and found to have an insignificant effect on the range of formation times, changing the original formation time range by less than the interpolation error of 0.01 Ma, but lowered the minimum radius size by slightly more than a kilometer (Table 1). As discussed in "Discussion", it is lower bound of thermal diffusivity's uncertainty that is responsible for the lowered minimum radius, so it Fig. 2 Effects of varying specific heat capacity. a The effect of uncertainty in specific heat capacity on the range of possible asteroid radii and formation times for an Itokawa-like parent body. Each of the three curves shown is actually made up of three isothermal lines corresponding to Itokawa's thermal constraints of having a peak temperature between 1073 and 1273 K and a central temperature greater than 973 K at 7.6 Ma, as in Figure 1. The solid black curve shows the original result, seen in Fig. 1. The dashed and dotted curves show the result when the specific heat capacity is increased or decreased, respectively. b, c Central temperature over time for two simulated bodies with end-member specific heat capacity values and formation time and radius values corresponding to the diamond markers on a. b A 21-km radius body formed at 1.78 Ma is able to narrowly avoid exceeding the maximum peak temperature constraint due to its greater heat capacity at 1.1 times the base value. c A 50-km body formed at 2.34 Ma scarcely fulfills the minimum peak temperature constraint due to its lower heat capacity at 0.9 times the base value can be said that its −10 % uncertainty resulted in a lowering of the minimum radius by more than 5%. Thermal diffusivity and specific heat capacity were also varied together, as their effects on the simulated bodies and on the overall results are not necessarily linear. Because there are four combinations of the end-member values of the uncertainties in these two parameters, four sets of simulations were performed and interpolated to create four sets of isothermal lines. The area of valid formation time and radius values is then found for each of these and combined to determine the total range of values allowed by this variation. This combined area is shown in Fig. 3. The variation of both parameters results in a formation time range identical to that found by varying only specific heat capacity, but the minimum radius value was found to be 19.0 km. This 1.5-km reduction in minimum radius is greater than that produced by the variation of either parameter separately (Table 1), but less than what would be expected if the reductions in minimum radius from varying the two parameters separately were simply added together.
The uncertainty in aluminum abundance could be explored the same way as thermal diffusivity and heat capacity, but it is not actually necessary to perform additional simulations to determine the effects of uncertainty in aluminum abundance, or in the magnitude of heat production in general, for this model. If the heat magnitude coefficient, A, of the heat production equation (Eq. 5) is changed by a multiplicative constant α , then this can be shown to be equivalent to a shift in time by an amount of time −τ ln(α) . Because the heat production equation is the only explicit function of time within our model, the thermal evolution of a body with a modified heat production magnitude is therefore equivalent to an identical body with unmodified heat production but a shifted formation time. The heat magnitude coefficient depends linearly on the wt.% of aluminum in the body material, and so the thermal evolution of a body with a varied aluminum abundance can be calculated by shifting the formation time of the thermal evolution of an identical body with the base aluminum abundance by −τ ln(α) , where α is the multiplicative variation from the base aluminum abundance. These bodies with a shifted formation time but unmodified aluminum abundance must still be simulated, but the creation of the isotherms over the free parameter space already requires numerous simulations with various formation times, so this does not incur additional computational cost. We used this shift in times to get the thermal evolution of bodies with varied aluminum abundances from the simulations already performed, and then created the combined isotherms and found the range of radius and formation time values in the manner already described.
The effects of varying only aluminum abundance proved significant for all of the uncertainty ranges considered. The smaller, conservative uncertainty range based on Sugiura and Fujiya (2014) resulted in an expanded formation time range of 1.72 to 2.40 Ma, but only lowered the minimum radius by less than a kilometer. That is, the ∼ ± 15% uncertainty in aluminum abundance resulted in about ± 7.5% changes in the bounds of the formation time range, which is nearly equivalent to the change in the formation time range per amount of uncertainty observed from specific heat capacity. The 2.5% reduction in minimum radius from the + 16% uncertainty in aluminum is notably less significant than the effects of heat capacity or thermal diffusivity on the minimum radius, per amount of uncertainty. The larger uncertainty range based on Jarosewich (1990) resulted in a large formation time range of 1.74-3.44 Ma, and a significantly lowered minimum radius size of 17.5 km. Here, the upper formation time bound is increased by more than 50% and the minimum radius is reduced by 15% due to the + 216% change in aluminum abundance. These relations between the uncertainty and the effects on the results are different from the proportions observed in the smaller uncertainty range due to the nonlinear nature of the effects of varying aluminum abundance. Lastly, varying only the aluminum abundance according to the range of Nakamura et al. (2014) resulted in a significantly earlier formation time range of 1.53-2.09 Ma, and an increased minimum radius of 21.1 km.
The effects of varying specific heat capacity, thermal diffusivity, and aluminum abundance together are then Fig. 3 Effects of varying multiple material properties. The range of possible radius and formation time values for Itokawa's parent body taking into account uncertainty in conductivity, heat capacity, and aluminum abundance. The properties are varied by the ranges of multiplicative factors given in Table 1. The variation of heat capacity and aluminum abundance lead to significantly larger ranges of possible formation times, with notably smaller minimum radius values, than when no variation or uncertainty in these parameters is considered Hallstrom and Bose Earth, Planets and Space (2023) 75:11 found by constructing the combined isothermal lines for each of the eight combinations of the end-member values of the parameters' uncertainty ranges and creating a single curve outlining the area of all of these ranges, as previously described ( Fig. 3; Table 1). When varying all three parameters, the overall minimum possible radius comes from the intersection of the 1273 K peak central temperature isotherm and the 973 K at 7.6 Ma central temperature isotherm generated using the highest aluminum abundance, lowest heat capacity, and lowest thermal diffusivity values considered. The latest possible formation time is found from the 1073 K peak central temperature isotherm from the same set of data, and the earliest formation time from the intersection of the 1273 K isotherm and the 973 K at 7.6 Ma isotherm generated using the lowest aluminum abundance, highest heat capacity, and highest thermal diffusivity values. Varying specific heat capacity and thermal diffusivity while varying aluminum abundance by the more conservative uncertainty range from Sugiura and Fujiya (2014) results in a total possible formation time range of 1.62 to 2.50 Ma after CAIs and a minimum radius size of 18.5 km. Using the larger but more tenuous uncertainty range from Jarosewich (1990) gives a much larger formation time range of 1.64-3.54 Ma with a minimum radius of 16.5 km. Varying specific heat capacity and thermal diffusivity along with the aluminum abundance range from Nakamura et al. (2014) gives a formation time range of 1.43-2.20 Ma with a minimum radius of 19.6 km.

Discussion
Varying aluminum abundance had a significant effect on the results ( Fig. 3; Table 1), even for our more conservatively created uncertainty range. For the uncertainty ranges considered here, aluminum abundance expanded the formation time range and lowered the minimum radius more than either of the other two parameters. However, this is primarily due to the relatively large amount of uncertainty determined for aluminum abundance in this case. Looking at the effects of varying each parameter per amount of variation or uncertainty, as done in Section 3, it can be seen that varying aluminum abundance has an equally significant impact on the formation time range as specific heat capacity and less significant impact on the minimum radius than thermal diffusivity. Aluminum abundance directly affects the total thermal energy available to a body, so reducing its abundance allows a body to form earlier without exceeding the maximum peak central temperature constraint, and increasing the aluminum abundance allows a body to form later while still reaching the minimum peak central temperature constraint. Additionally, when a body is formed later it has less time to cool before 7.6 Ma, meaning bodies with smaller radii can retain enough heat to meet the minimum temperature constraint at that time.
The variation of specific heat capacity also significantly impacted our results (Fig. 2), despite its smaller range of variation, in agreement with previous studies on the importance of specific heat capacity in thermal evolution models, notably Ghosh and McSween Jr. (1999). The specific heat capacity determines the magnitude of the temperatures that can be reached by the material for a given amount of thermal energy, and so affects the results in a similar fashion as varying the aluminum abundance. An increased specific heat capacity allows a valid Itokawa parent body to form earlier, and a decreased heat capacity allows such a body to form later and smaller.
While our defined uncertainty range for thermal diffusivity was relatively large, the variation of the parameter had no noticeable effect on the formation time range. The −10 % uncertainty in thermal diffusivity did decrease the minimum radius by more than 5%, however, which is a more significant reduction in minimum radius per amount of uncertainty than observed by varying aluminum abundance or specific heat capacity. For the types of model parameters and bodies considered here, the rate at which heat is being generated within the material is much greater than the rate at which heat can be conducted to the surface during the initial heating phase. The peak temperatures reached in the centers of our simulated bodies are therefore not strongly dependent on the thermal conductivity or diffusivity of the material. However, the 7.6-Ma minimum temperature constraint is noticeably affected by thermal diffusivity, as the lower value of thermal diffusivity allowed smaller bodies to be valid parent bodies of Itokawa when they would otherwise have failed to maintain the minimum temperature needed at 7.6 Ma, which is shown by the lowered minimum radius of the varied results. It is also of note that while varying specific heat capacity and thermal diffusivity both lead to reductions in the minimum radius, the amount by which the minimum radius is lowered when the two parameters are varied simultaneously, 1.5 km, is slightly less than the sum of the reductions observed when the parameters are varied separately, 0.5 km and 1.2 km (Table 1).
In addition to these three material parameters, the effect of uncertainty in the initial and surface temperatures was briefly investigated. We performed several simulations of a 30-km radius body formed at 2.2 Ma where initial and surface temperature were varied from 150 to 300 K and found that the peak central temperatures were all within 100 K of the original value, and so these parameters were not investigated further.
To summarize, the parent body of Itokawa would have likely accreted between 1.6 and 2.5 Ma after time 0, but formation times as early as 1.4 Ma or as late as 3.5 Ma are possible if outlier aluminum abundance measurements from Nakamura et al. (2014) or Jarosewich (1990), respectively, are considered. The range of allowed formation times from our model is larger than previously estimated for Itokawa or OCs generally. Previous estimations include the 1.9-2.2 Ma range for a > 20 km radius Itokawa-like parent body found by Wakita et al. (2014), the 1.6-1.8 Ma range for a 20-40 km radius L chondritelike parent body from Doyle et al. (2015), or 2.04-2.24 Ma for a ∼ 50-km radius ordinary chondrite-like parent body from Sugiura and Fujiya (2014). The optimized formation times between 1.8 and 1.9 Ma for the > 100 km radius H and L chondrite parent bodies by Henke et al. (2012b) and Gail and Trieloff (2019), respectively, are also both within our determined range. If the formation of OCs is very late (3.5 Ma; Table 1), then it implies rapid thermal metamorphism of OC parent bodies and this could be a case for Itokawa-like bodies. Additionally, the possibility that the radius of the Itokawa parent body could indeed be smaller than 20 km exists according to the bodies allowed by the criteria used here, specifically due to the effects of reduced thermal diffusivity ( ∼ 10% lower), greater aluminum abundance (> 16% higher), or both, contrary to earlier results.

Importance of measuring physical properties of meteorites
The results shown here stress the importance of accurate data of the material properties of chondritic meteorites. The database of heat capacity of OCs is mostly limited to H and L and not LLs, and the majority of measurements on the thermal properties of chondritic material is limited to temperatures of 500 K or lower despite a large number of chondrites experiencing larger temperatures, greater than 1000 K in most instances, e.g., in the centers of the H (Henke et al. 2012b) and L (Gail and Trieloff 2019) parent bodies. The thermal properties of OCs (and all other meteorite types) are functions of their composition and mineralogy, and so additional measurements of a wide range of OCs are needed for thermal evolution models such as these to be sufficiently accurate. There are added complications in using meteorite samples to understand asteroid populations, as measurements of the material may be affected by the long residence times of the meteorite on earth, fragile components of the material could have been destroyed or altered over time, and the material could display different behavior in earth's environment than in the environment of an early solar system parent body. In spite of these difficulties, such measurements of physical properties are needed, especially now for all types of carbonaceous chondrite meteorites within the CM and CR classes owing to the Ryugu samples in our collection, and the Bennu samples that are scheduled to arrive soon.
There are sources of error in our model that remain to be quantitatively and thoroughly analyzed. Some of these include the omission of additional heat sources (e.g., additional radionuclides like 60 Fe or 40 K), the simplifying assumptions that the body is perfectly spherical and is instantaneously formed, and the omission of physical phenomena such as those arising from material porosity and heterogeneity (as discussed in more complex models such as that of Henke et al. (2012a)). Both the uncertainties in material parameters and these additional sources of error exist to some degree in almost all small-body thermal evolution models. While the relative magnitude of these errors depends on the specific construction and purposes of each such model, the significance of uncertainty in the material properties discussed here is likely to have similar behavior and importance in other smallbody models.

Conclusions
(1) Prior to this work thermal evolution models and previously published thermometry data limited the accretion time of L and LL OCs to between 1.9 and 2.2 Ma (Wakita et al. 2014), 2.0-2.2 Ma (Sugiura and Fujiya 2014), 2.05-2.25 Ma (Blackburn et al. 2017), or 1.8-1.9 Ma (Gail and Trieloff 2019). This study uses uncertainties in key bulk properties of OC materials to allow an Itokawa-like parent body to have accreted anytime between 1.6 and 2.5 Ma after the formation of CAIs with a radius less than 20 km. Furthermore, if outlier measurements of aluminum abundance in Itokawa particles  or LL chondrites (St. Mesmin in Jarosewich (1990)) are considered, the formation time could be as early as 1.4 Ma or as late as 3.5 Ma. Uncertainties in aluminum abundance and specific heat capacity both contributed to the new uncertainty in the formation time of the Itokawa-like parent body, by similar amounts of time per amount of uncertainty. Uncertainty in thermal diffusivity primarily resulted in uncertainty in the radius of the parent body, to a greater extent than either aluminum abundance or specific heat capacity.
(2) Accurate and precise laboratory measurements of chondritic material properties would greatly improve the efficacy of these models, and measurements of temperature-dependent properties such as specific heat capacity and thermal diffusivity at temperatures up to the peak temperatures reached in the body of that particular material are required for accurate thermal modeling of undifferentiated asteroid parent bodies.
(3) Interpolating the results of simulated bodies can be used to characterize the formation time and radius values which comply with known thermal expectations and constraints. For computationally expensive simulations, which many thermal evolution models are, trying to account for uncertainty in multiple parameters through model variations would likely require an interpolation approach to be reasonably pursued. We obtained sufficiently smooth and consistent isothermal lines over the space of the free parameters of radius and formation time using approximately 80-120 simulated bodies over the space of interest for each combination of varied parameters investigated. As here, care should also be taken to ensure that it is only the end-member variations of uncertain properties that produces results of interest, which also allows for a reduced number of required simulations compared to a more thorough use of the entire variation or uncertainty ranges.

OC Ordinary chondrite CAI Calcium aluminum inclusion Ma
Million years