Dust environment of active asteroids P/2019 A4 (PANSTARRS) and P/2021 A5 (PANSTARRS)

We report on the characterisation of the dust activity and dynamical evolution of two faint active asteroids, P/2019 A4, and P/2021 A5, observed with the 10.4m GTC using both imaging and spectroscopy. Asteroid P/2019 A4 activity is found to be linked to an impulsive event occurring some $\pm$10 days around perihelion, probably due to a collision or a rotational disruption. Its orbit is stable over 100 Myr timescales. Dust tail models reveal a short-term burst producing (2.0$\pm$0.7)$\times$10$^6$ kg of dust for maximum particle radius rmax=1 cm. The spectrum of P/2019 A4 is featureless, and slightly redder than the Sun. P/2021 A5 was active $\sim$50 days after perihelion, lasting $\sim$5 to $\sim$60 days, and ejecting (8$\pm$2)$\times$10$^6$ kg of dust for rmax=1 cm. The orbital simulations show that a few percent of dynamical clones of P/2021 A5 are unstable on 20-50 Myr timescales. Thus, P/2021 A5 might be an implanted object from the JFC region or beyond. These facts point to water ice sublimation as the activation mechanism. This object also displays a featureless spectrum, but slightly bluer than the Sun. Nuclei sizes are estimated in the few hundred meters range for both asteroids. Particle ejection speeds ($\sim$0.2 m/s) are consistent with escape speeds from those small-sized objects.


INTRODUCTION
Active asteroids constitute a new class of objects in the solar system. They are characterised by being located in the main asteroid belt, but, contrary to most objects in the belt, display cometary appearance, i.e., dust comae and tails. Some 40 objects of this kind have been discovered so far. The reasons of their activity are rather diverse, including impactinduced, rotational break-up, thermal fracture, or ice sublimation (Jewitt, Hsieh, & Agarwal 2015). While most of those objects seem native to the main belt, some of them have been shown to become unstable on timescales of a few tens of Myr (e.g. Haghighipour 2009;Hsieh & Haghighipour 2016). Given the variety of phenomena that might lead to dust ejection, it is very convenient to increase the sample statistics to provide a better knowledge of the physics involved. To that end, current sky surveys such as the Panoramic Survey Telescope and Rapid Response System (PANSTARRS) and facilities to come such as the Vera Rubin Observatory are of the utmost importance in detecting these typically faint objects.
Asteroid P/2019 A4 was discovered on UT 2019 Jan. 10.4 by PANSTARRS ). Its orbital elements a=2.614 au, e=0.0896, and i=13.32 • yield a Tisserand pa-⋆ E-mail: fernando@iaa.es rameter respect to Jupiter of T J =3.36, as most asteroids, yet displaying a sizeable coma and tail. It is located in the middle portion of the belt, as P/2010 A2 (e.g., Jewitt et al. 2010), P/2016 G1 (e.g., Moreno et al. 2016;Hainaut et al. 2019), or (6478) Gault (e.g., Moreno et al. 2019;Jewitt et al. 2019), but with a remarkably small eccentricity. In these three cases, the activity has been found to be linked to short-term events, either a rotational mass loss or an impact as the possible causes of the events. P/2021 A5 was also discovered by PANSTARRS (Weryk 2021), on UT 2021 Jan. 06, when displaying a condensed coma and a 4 ′′ tail. Its orbital elements (a=3.047 au, e=0.14, i=18.19 • ) give T J =3.147, and place the object close to the outer main belt where ice sublimation presumably starts to become dominant over other mechanisms. Clear examples of this mechanism at play are 324P/La Sagra (e.g., Moreno et al. 2011;Hsieh & Sheppard 2015;Jewitt et al. 2016) and 133P/Elst-Pizarro (Hsieh, Jewitt, & Fernández 2004;Hsieh et al. 2010), that have shown recurrent activity.
In this work, we present images and spectra of P/2019 A4 and P/2021 A5 obtained with the Gran Telescopio Canarias (10.4m aperture) on the island of La Palma. We study the dynamical evolution of the objects to shed some light on their origin and on the causes of their activity, and apply Monte Carlo dust tail models to the images obtained to determine the dust properties, and to place further constraints on the activation mechanisms.

ORBITAL DYNAMICS SIMULATIONS
To asses whether the two objects under study are native to the main belt or interlopers coming from elsewhere, we propagated their orbits backward in time up to 100 Myr. We remark that the short arc of observations of these objects (in particular, P/2021 A5) results in relatively large uncertainties in their best-fitted orbital parameters. In any case, for each object, we integrated the orbits of 200 dynamical clones according to the statistical uncertainty of the current orbital elements. The orbital elements of those clones were generated using the covariance matrix (e.g. Milani and Gronchi 2010). The six-component vector of orbital elements x ′ of the dynamical clones are calculated according to the expression where x is the six-component vector of the nominal orbital elements (best-fit solution), and A is a matrix verifying AA ⊺ = C, where C is the covariance matrix. This matrix is obtained from the JPL Small-Body Database for each asteroid, and matrix A is obtained from C by a Cholesky decomposition using the FORTRAN implementation as described in Press et al. (1992). Ψ is a six-dimensional vector whose components are normally distributed (Gaussian) deviates with zero mean and unit variance, which are also obtained from the corresponding routine described in Press et al. (1992).
The time integrations were performed using the Bulirsch-Stöer integrator of MERCURY package (Chambers 1999). All the eight planets were included as major bodies, while all the dynamical clones were considered as massless particles. Nongravitational forces were not included. We used an integration time step of 10 days. The orbits of all 200 dynamical clones generated for P/2019 A4 were stable over the 100 Myr integration time, meaning that this object is very likely native from the main belt. The case of P/2021 A5 is very different, however, because of its proximity to the 9:4 Jupiter resonance (a=3.031 au). We have noticed that 8 out of the 200 P/2021 A5 clones become eventually unstable, being either ejected from the solar system or experience a collision with some of the inner planets or the Sun on timescales of order 20-50 Myr. Some of those ejected clones have intermediate JFClike and even Centaur-like orbits before ejection. Figure 1 displays intermediate orbital elements and the evolution of the Tisserand parameter of some of those unstable clones at 1000 year-intervals. The uppermost panel shows the typical dynamical evolution of the unstable clones, where the 9:4 resonance induces excitation in eccentricity and inclination, driving the particle to Mars and Jupiter crossing orbits, and spending some time in JFC-like orbits before being ejected from the solar system. Among those unstable clones, we have found two atypical cases in the dynamical evolution. In one of those (central panels in Figure 1), the particle migrates inwards, experiencing several episodes of resonant oscillations, and being finally trapped in the strong 5:2 Jupiter resonance region at 2.825 au, from where its eccentricity is increased up to almost unity until a collision in the inner solar system occurs. The other peculiar case (lowermost panels in Figure 1) correspond to a particle which, after being excited by the 9:4 resonance, spends some time in the JCF region, and then experiences several resonant episodes in the Centaur region beyond Jupiter's orbit before being finally ejected from the solar system. While these calculations are not intended to provide robust statistical results, as the number of test particles is not sufficient, and the nominal orbital parameters have moderately large uncertainties, they do indicate possible different origins of these two active asteroids under analysis. As shown by the dynamical experiments, P/2019 A4 is very likely a native member of the main belt, while for P/2021 A5 there is some probability that it comes from elsewhere, possibly from the JFC or Centaur regions, as a result from the dynamical evolution of the clone orbits. This suggests that P/2021 A5 might be an ice-bearing object. In this regard, we note that possible dynamical pathways from the JFC region to the main belt have been previously explored by Fernández, Gallardo, & Brunini (2002), who did not find any, and more recently by Hsieh and Haghighipour (see Hsieh & Haghighipour 2016, and references therein), who concluded that the number of JFC-like interlopers in the main belt, albeit likely small, might be non-zero.

OBSERVATIONS
Images and spectra of asteroid P/2019 A4 have been obtained under photometric conditions on 2019 February and 2019 March, and of P/2021 A5 on 2021 February. The log of observations is shown in Table 1. Images of P/2019 A4 were obtained on a CCD camera using a Sloan r ′ filter with the Optical System for Image and Low Resolution Integrated Spectroscopy (OSIRIS) camera-spectrograph (Cepa 2010) attached to the 10.4m Gran Telescopio Canarias (GTC) at the Roque de los Muchachos Observatory on the island of La Palma (Spain). For P/2021 A5, additional images using the same setup were also obtained through g ′ and i ′ filters. The image scale was 0.254 ′′ pixel −1 . The images were biassubtracted, flat-fielded, and flux calibrated using standard procedures. Median stack images were computed from the available frames on each night. The set of images of these faint asteroids at the r ′ filter are displayed in Figure 2. The seeing was near 1.25 ′′ FWHM during those nights, i.e., about 5 pixels FWHM.
The spectral images were bias and flat-field corrected, using lamp flats. Spectra were extracted and collapsed to one dimension, using an aperture of ±5 pixels centred at the maximum of the intensity profile of the object. The wavelength calibration was performed using Xe+Ne+HgAr lamps. Finally, the 3 spectra obtained of each comet were averaged to obtain the final spectrum. To obtain the normalized reflectance spectra of P/2019 A4 and P/2021 A5 shown in Figure 3 we observed two G2V stars that were used as solar analogues. For P/2019 A4, we used the standard SA 98 978 (Landolt 2009) , while for P/2021 A5, we used SA 93 101 (Galadí-Enríquez, Trullols, & Jordi 2000). The colours of those stars along with the solar colours from Ramírez et al. (2012) are given in Table 3.
The reflectance spectra are featureless, with no evidence of CN, C 2 or C 3 emission bands typical of coma gas emission nor the wide absorption band typical of the stony class Table 1. Log of the observations. R and ∆ denote the asteroid heliocentric and geocentric distances, respectively. χ is the plane angle, i.e., the angle between the observer and the asteroid orbital plane, and α is the phase angle.  Table 2. Photometric results. Sloan magnitudes obtained using a 10 pixels diameter apertures are given together with the lower limit of the absolute magnitude and the derived colours. asteroids. Unfortunately, the SNR of both spectra is too low to attempt to derive reliable gas production rates. The P/2019 A4 spectrum is slightly red with slope S ′ = 4±1%/100 nm typical of a primitive X-type asteroid. The P/2021 A5 spectrum is very noisy, and slightly blue with a slope S ′ = −3 ± 1%/100 nm as typical of a primitive B-type asteroid. From the computed r ′ , g ′ , and i ′ magnitudes (see Table 2), we get (g ′ −r ′ )=0.63±0.05 and (r ′ −i ′ )=0.06±0.05. The colours of the Sun are (g ′ − r ′ ) ⊙ =0.44 and (r ′ − i ′ ) ⊙ =0.11. The measured colour index (r ′ − i ′ ) is therefore consistent with the slightly blue colour respect to the Sun observed in the spectrum. In addition, we have obtained the colour index (r ′ − i ′ ) from the P/2021 A5 binned spectrum, resulting in (r ′ − i ′ )=0.09, which is within the error bars of the measured index from image photometry.
For modelling purposes, we take the geometric albedo of these objects as p v =0.106 for P/2019 A4 (see Mainzer et al. 2012, median value for X-type asteroids), and p v =0.07 for P/2021 A5 (see Alí-Lagoa et al. 2016, averaged value for Btypes). The density of these bodies were taken as 1850 kg m −3 and 2380 kg m −3 for P/2019 A4 and P/2021 A5, respectively (Carry 2012). The corresponding linear phase coefficients are calculated from the magnitude-phase relationship by Shevchenko (1997) as 0.036 mag deg −1 and 0.04 mag deg −1 for P/2019 A4 and P/2021 A5, respectively.
In order to set a stringent upper limit to the asteroid radii, we converted the absolute magnitudes observed to diameters using the usual relationship from Pravec & Harris (2007). We remark that these diameters are just upper limits that we will use later in the mode-ling procedure, given the strong contamination of dust around those objects. Absolute magnitudes in the r ′ filter are computed as H ′ r = r ′ − 5 log(R∆) − ηα, where η is the linear phase coefficient. Assuming a (B −V ) colour index approximately solar (see Ta-

DUST MODELLING
We proceed in the same way as in previous works, were we made use of our Monte Carlo dust tail code to characterise the dust environment of comets and active asteroids. Recent applications of the model can be found in Moreno et al. (2019) and de León et al. (2020). The idea is to build up a dust tail for a certain object at a given time given a series of input parameters such as the size distribution of the particles, dust mass loss rate, and ejection speeds. The particles are assumed to be characterised by their density, the geometric albedo, and the linear phase coefficient. In the Monte Carlo procedure, the tail brightness is computed by adding up the contribution to the brightness of each particle ejected from the comet or asteroid nucleus. The trajectory of each particle depends on the ratio of radiation pressure force to the gravity force, i.e., the β parameter, which is defined as β = C pr Q pr /(2ρr), where C pr =1.19× 10 −3 kg m −2 is the radiation pressure coefficient, Q pr is the scattering efficiency for radiation pressure, and ρ is the particle density. Q pr is taken as 1, as it converges to that value for absorbing particles of radius r 1 µm (see e.g. Moreno et al. 2012, their Figure 5). Since the number of physical parameters is large, several assumptions must be made to make the problem tractable. We assume the geometric albedo, linear phase coefficient, and density, as estimated in the previous section for each asteroid. In addition, the particles are assumed to be distributed in a broad power-law size distribution function n(r) ∝ r κ , where κ is the power-law index.
If the dust production rate is not sufficiently high as to dominate the scattering cross section, the nucleus surface might have a significant contribution to the observed brightness. For the two targets under study, this might be the case, so we have added up the nucleus contribution to the brightness by considering its effective cross section and the geometric albedo and linear phase coefficient given in the previous section. The nuclear radii (R N ) are then considered as one of the free parameters of the model, always subjected to the constraint of being smaller than the upper limits estimated in Section 2 from the absolute magnitudes.
The remaining dust parameters involved in the dust tail brightness computation are taken as free parameters to be fitted by a multidimensional fitting algorithm, namely the downhill simplex method of Nelder & Mead (1965), which has been implemented in FORTRAN language by Press et al. (1992). The best-fit parameters are found by minimising the squared sum of the differences between the modelled and measured tail brightness for the GTC images. To perform an appropriate comparison between observed and modelled tails, each modelled tail is convoluted with a Gaussian function having a FWHM equal to the measured seeing. The fitting parameters are the power-law exponent of the size distribution κ, the ejection speeds, and the dust mass loss rate. To keep the number of free parameters to a minimum, the dust mass loss rate function is assumed as a Gaussian function with parameters M t (the total dust mass ejected), the time of maximum dust loss rate (t 0 ) and the FWHM of the Gaussian, which gives a measure of the effective time span of the emission event. The peak dust loss rate, (dM/dt) 0 , is related to the total dust mass loss and the FWHM through the equation M t = 1.06(dM/dt) 0 FWHM. The ejection speed is assumed to follow the equation v = v 0 β γ , where v 0 and γ are fitting parameters of the model. Then, the set of seven fitting parameters (NP = 7) is M t , t 0 , FWHM, v 0 , γ, κ, and the nucleus radius, R N . To start the execution of the code, an initial simplex must be set with NP+1=8 sets of input parameters that we choose to vary broadly between reasonably expected minimum and maximum limits. Since the best-fit set of parameters found in the downhill simplex method necessarily corresponds to a local minimum of the fitting function, we repeated the procedure for a variety of input starting simplex parameters in an attempt to find the lowest of those local minima.
Some information on the ejection event timings and the size of the ejected particles can be retrieved from the syndynesynchrone network associated to each observation date. For asteroid P/2019 A4, we show the syndyne-synchrone map associated to the image taken on 2019 March (see Figure 4). This network suggests that dust emission must have occurred near perihelion, and likely concentrated ∼ ±20 days respect to that date. In addition, particles smaller than about 30 µm do not contribute significantly to the tail brightness, although one must always keep in mind that the syndyne-synchrone map refer to the geometric loci of particles ejected with zero velocity from the nucleus. For P/2021 A5, the syndynesynchrone map for the observation date is plotted in Figure 5. In this case, we see that the activity must have occurred when the object was approaching perihelion and later. Owing to the geometry of the observation from Earth, post-perihelion synchrones are closely spaced becoming more difficult in this case  to predict an event date from the syndyne-synchrone network. The syndyne curves indicate that dust particles smaller than about 30 µm contribute very little to the observed tail.
We use these syndyne-synchrone networks to establish safe limits in the starting simplex for the time interval during which the asteroids were active, in terms of t 0 and FWHM of the Gaussian defining the dust loss rate profile, as well as the minimum and maximum particle radii of the size distribution function, which we set to 10 µm and 1 cm, respectively, for both asteroids.

RESULTS AND DISCUSSION FOR P/2019 A4
For asteroid P/2019 A4 we found the best-fit parameters as shown in Table 4, and the observed and simulated tails in  Figure 6. The FWHM is always constrained to FWHM<20 days, and the total dust mass loss is given by M t =(2.0±0.7)×10 6 kg. This short duration event indicates that the most probable cause is either an impact or a rotational disruption. Although the time of peak activity is very close to asteroid perihelion, we rule out, in principle, ice sublimation as a possible mechanism because of the asteroid position in the central part of the belt, a too hot location for stable reservoirs of water ice to exist, and its orbital stability over a timescale of 100 Myr or longer. In addition, the low eccentricity of this object makes ice sublimation at perihelion even less likely.
The nominal ejection speeds (v = 1.6β 0.2 m s −1 ) are found to depend slightly on the particle radius (the exponent, γ, is γ=0.2), ranging from a minimum of 0.2 m s −1 to a maximum of 0.8 m s −1 for particles of radii 1 cm and 10 µm, respectively. The derived nominal nucleus radius, R N =170 m, implies an escape speed of v esc =0.23 m s −1 which is very consistent with the ejection speed of the largest particles considered (r=1 cm).
If an impact were the cause of the observed activity, the impactor's size can be roughly estimated from the ejected mass following an argument similar to that by Jewitt, Ishiguro, & Agarwal (2013) for asteroid P/2010 A2. For an average collision speed of ∼5 km s −1 in the main belt, and an ejecta escape velocity of 0.23 m s −1 , the ratio of ejecta mass to projectile mass, assuming that impactor and target have the same density, is of order M e /M p ∼10 4 (Housen & Holsapple 2011). Then, if the ejected mass is M e = 2 × 10 6 kg, the impactor mass becomes M p =200 kg, corresponding to a spherical object of only 0.3 m in radius (ρ=1850 kg m −3 ). On the other hand, a rotational disruption might have also occurred: the low speed of the particles ejected and the fact that the event might last as much as 20 days point to that possibility. Besides, the timescale for Yarkovsky-O'Keefe-Radzievskii-Paddack (YORP) induced rotational acceleration for objects having radius 6 km is always shorter than the collision timescale (Jacobson et al. 2014). Thus, for an asteroid at 2.5 au from the Sun and a YORP coefficient of Y =0.01 (which is related to the asteroid shape), τ YORP = 42R 2 A Myr, where R A is the radius of the asteroid in km (see Jacobson et al. 2014), resulting in τ YORP =1.2 Myr for P/2019 A4. The collisional timescale for the nominal asteroid radius (170 m) is τ coll ∼140 Myr (see Bottke et al. 2005, their Fig. 14), i.e., more than two orders of magnitude longer. However, at this point we are unable to favour one of these two mechanisms. Observations of the lightcurve of the target could help in solving the problem. Nevertheless, the faintness of the object might preclude any future attempt to measure its rotation period.
6 RESULTS AND DISCUSSION FOR P/2021 A5 Table 4 gives the results for the best-fit parameters concerning asteroid P/2021 A5, and the observed and modelled im-ages in Figure 7. In this case, the event time is shifted by +48 days (nominally) after perihelion passage, but with a large uncertainty of ±30 days. The event duration is constrained between 5 and 60 days, and the total dust mass loss released is M t =(8±2)×10 6 kg. The nominal particle ejection speeds (0.15 m s −1 ) are found to be independent of size (γ=0.0). The escape speed corresponding to the derived R N =150 m nucleus is 0.17 m s −1 , which is of the same order of the ejection speeds. Regarding the activation mechanism(s), the orbital dynamics of a few percent of the dynamical clones point to a JFC origin for this object and, consequently, it might be an ice-bearing asteroid. The ice-driven activity is compatible with the duration of the activity, FWHM=39 days (nominal). However, mass shedding from rotational instabilities is also feasible. Perhaps the less probable cause is an impact, because of the lower limit of the duration of FWHM>5 days, although it cannot be ruled out. Further observations of this target in future apparitions could shed light on the responsible mechanism(s).

CONCLUSIONS
Key conclusions on the observation and modelling of active asteroids P/2019 A4 and P/2021 A5 are as follows.
1) Orbital dynamics simulations show that P/2019 A4 moves in a stable orbit over long (100 Myr) timescales. Asteroid P/2021 A5 is located close to the 9:4 Jupiter resonance region and the dynamical evolution of 4% of its clones indicates a possible JFC origin.
2) The spectrum of P/2019 A4 is slightly red with slope S ′ = 4 ± 1%/100 nm typical of a X-type asteroid while the of P/2021 A5 spectrum is slightly blue with a slope S ′ = −3 ± 1%/100 nm, typical of a primitive asteroid of B type.
3) From the Monte Carlo dust tail modelling of GTC images on these targets, we estimate that P/2019 A4 was active for a maximum period of 20 days (FWHM), which indicates, probably, an activation mechanism related to either a rotational disruption or an impact. This is also supported by the stability of its orbit in the central part of the belt. This increases the statistics of bodies in the middle belt whose activity is reported as short-lived. Regarding P/2021 A5, the longer possible activity period up to 60 days, and the possibility of a JFC origin as revealed by the orbital dynamics simulations could indicate an ice-driven activity, but we are unable to rule out other hypotheses due to the large uncertainty (plausible FWHM values in the 5-60 days range) in the estimated duration.
4) The total dust mass ejected from these two asteroids are, for a maximum particle radius ejected of 1 cm, (2.0±0.7)×10 6 kg and (8±2)×10 6 kg for P/2019 A4 and P/2021 A5, respectively, and the ejection is concentrated close to the perihelion passage for P/2019 A4, and possibly shortly after perihelion passage for P/2021 A5, although with an uncertainty of ±30 days.
5) The derived ejection speeds for both targets ∼0.2 m s −1 are consistent with the escape speeds of the nuclear radii estimated from the Monte Carlo modelling, which turned out to be in the range of 100-240 m (nominal value 170 m) for P/2019 A4, and smaller than 500 m (nominal value 150 m) for P/2021 A5.  Table 4. Upper panels: observation (a) and model (b) images on 2019 February. The rightmost panel (c) displays the observed (black contours) and modelled (red contours) isophotes. Innermost isophote corresponds to 2×10 −14 solar disk intensity units. Isophotes decrease in factors of two outwards. Lower panels: same as the upper panels but for the 2019 March image. Isophotes decrease in factors of two outwards, the innermost one corresponds to 5×10 −15 solar disk intensity units.

DATA AVAILABILITY
The spectral and image data shown in this manuscript are available from the GTC Public Archive after the one-year proprietary period is over. Figure 7. Observation and best fit model for P/2021 A5, as indicated. The best-fit model parameters are the nominal parameters in Table 4. In the lowermost panel, the observed (black contours) and modelled (red contours) isophotes are displayed. Innermost isophote corresponds to 3×10 −14 solar disk intensity units. Isophotes decrease in factors of two outwards.