The Climate and Compositional Variation of the Highly Eccentric Planet HD 80606 b -- the rise and fall of carbon monoxide and elemental sulfur

The gas giant HD 80606 b has a highly eccentric orbit (e $\sim$ 0.93). The variation due to the rapid shift of stellar irradiation provides a unique opportunity to probe the physical and chemical timescales and to study the interplay between climate dynamics and atmospheric chemistry. In this work, we present integrated models to study the atmospheric responses and the underlying physical and chemical mechanisms of HD 80606 b. We first run three-dimensional general circulation models (GCMs) to establish the atmospheric thermal and dynamical structures for different atmospheric metallicities and internal heat. Based on the GCM output, we then adopted a 1D time-dependent photochemical model to investigate the compositional variation along the eccentric orbit. The transition of the circulation patterns of HD 80606 b matched the dynamics regimes in previous works. Our photochemical models show that efficient vertical mixing leads to deep quench levels of the major carbon and nitrogen species and the quenching behavior does not change throughout the eccentric orbit. Instead, photolysis is the main driver of the time-dependent chemistry. While CH$_4$ dominates over CO through most of the orbits, a transient state of [CO]/[CH$_4$}] $>$ 1 after periastron is confirmed for all metallicity and internal heat cases. The upcoming JWST Cycle 1 GO program will be able to track this real-time CH$_4$--CO conversion and infer the chemical timescale. Furthermore, sulfur species initiated by sudden heating and photochemical forcing exhibit both short-term and long-term cycles, opening an interesting avenue for detecting sulfur on exoplanets.


INTRODUCTION
The orbital eccentricity of planetary systems sheds light on formation history and dynamical evolution. The demographic of eccentricity has been constructed through radial velocity (e.g., Shen & Turner 2008), transit (e.g., Shen & Turner 2008) surveys, and synthesis analysis (e.g., Kipping 2013). The majority of close-in planets have circularized orbits owing to the strong tidal interaction (Pont et al. 2011), which dissipates eccentricity on a timescale that steeply shortens with the semi-major axis at a power-law rate (Goldreich & Soter 1966). On the other hand, the discovered exoplanets exhibit a diverse eccentricity distribution beyond a semi-major axis of ∼ 0.1 AU (Butler et al. 2006;Kane et al. 2012). To understand the mechanisms causing high eccentricity and its impact is of great interest. Kane & von Braun (2009) have also shown that eccentric planets on average have increasing transit and eclipse probability as compared E-mail: shang-min.tsai@physics.ox.ac.uk to their circular-orbit counterparts with the same orbital period, making them ideal observation targets.
The Shoemaker-Levy 9 impact event in 1994 presented a dramatic but one-time example of what can temporal evolution teach us (Harrington et al. 2004;Hammel et al. 2010). Individually, an eccentric planet manifests climate variability periodically, as the stellar irradiation varies across the orbit. As opposed to the generally steady state of tidally-locked planets in circular orbits, the eccentricity induced variability provides unique information on how physical and chemical processes interact in the atmosphere. Furthermore, extremely eccentric planets could spend most of the orbit receiving Earth-like irradiation while being briefly heated near periastron and affording direct observations. This climate regime transition provides an illuminating connection between Solar System planets to close-in exoplanets.
HD 80606 b is a gas giant (4 M J ) in the binary system of HD 80606 and HD 80607. The extremely high eccentricity of HD 80606 b (e 0.933) places it among rare handful of planets with e greater than 0.8. The high eccentricity and close orbit of HD 80606 b are potentially a result of "Kozai migration", i.e. the companion star perturbed the planet with an initially inclined orbit into an eccentric orbit, then tidal dissipation drew the planet inward after Kozai oscillation stopped (Wu & Murray 2003). During the planet's close encounters, its planetary spin is expected to be quickly synchronized with the orbital revolution at periastron by the strong tidal force (pseudo-synchronous state). However, Spitzer observations have suggested that the exact rotational period might differ from the pseudo-synchronous state (de Wit et al. 2016;Lewis et al. 2017).
The eccentric configuration and the variation in stellar radiation of the HD 80606 system are illustrated in Figure 1. The unique dynamics of HD 80606 b have attracted various theoretical and observational studies. Laughlin et al. (2009) provided the first look at the Spitzer photometric observations at 8 µm. By measuring the increase of the planetary flux, Laughlin et al. (2009) estimated the global radiative timescale to be ∼ 4.5 hours, much shorter than that in cooler atmospheres in Solar System (∼ days). de Wit et al. (2016) followed up with the 4.5 µm channel observation spanning a longer phase coverage around the periastron. de Wit et al. (2016) further considered the effect of the revolving substellar longitude by planetary rotation and inferred a rotational period longer than the pseudo-synchronous period (∼ 93 hours). The tidal quality factor of HD 80606 b is also estimated to be much higher than that of Jupiter, indicating a lower rate of tidal dissipation that is consistent with the high eccentricity the system retained. Finally and excitingly, HD 80606 b will soon be observed by the NIRSpec and MIRI instruments on JWST during Cycle 1 General Observers (GO) program JWST- GO-2008 (PI Kataria) and JWST-GO-2488 (PI James Sikora). These observations will provide time-series spectra to reveal how the planet responds to the flash heat around periastron in detail.
From the modeling effort, 1D time-stepping radiative models have been applied to studying the thermal evolution of the full phase (Iro & Deming 2010;Mayorga et al. 2021). To probe the atmospheric dynamics, Langton & Laughlin (2008) applied a sallow-water model and identified the fast zonal flow driven by rapid heating near periastron. Lewis et al. (2017) explored the effects of rotational period with a 3D GCM and found that the timing of the planetary flux strongly depends on rotation. Lewis et al. (2017) also suggested that observed phase curve near periastron could be explained by the emergence and dissipation of optically thick clouds.
While most of the previous works about eccentric planets focused on the thermal response, they often assume a constant composition when analyzing the light curves. However, changes in the atmospheric chemistry can significantly alter the photospheric levels probed by the observations (Dobbs-Dixon & Cowan 2017;Parmentier et al. 2021) and hence the derived planetary properties (e.g., rotational rate). The assumption of maintaining thermochemical equilibrium instantaneously is also not expected to hold due to various disequilibrium processes, such as atmospheric mixing, photochemistry, and the stellar flux variation along the eccentric orbit, etc. Visscher (2012) provides an analytical overview of the chemical conversion on eccentric planets and indicates that the disequilibrium processes play an important role in the CO-CH 4 interconversion.
Despite the great effort devoted to study HD 80606 b, a comprehensive understanding of how chemistry and dynamics evolve in an eccentric system is still lacking. To unveil the essential atmospheric properties en masse, we assemble 3D GCM and photochemical models in this work to study the climate and compositional response in tandem. In Section 2, we discuss GCM modeling and the global thermal and dynamical responses to the rapid heating. In Section 3, we introduce the 1D time-dependent photochemical model and show the chemical response from simple timescale comparison T eq (K) Figure 1. The configuration of the HD 80606 system. The eclipse illustrates the orbit of HD 80606 b, color coded by the planet's equilibrium temperature. The black circles mark the secondary eclipse and primary transit and grey circles mark periastron and apoastron. The dotted circle around the star shows the circular orbit with the same rotational period as that in a pseudo synchronous rotation adopted in our model. to detailed photochemical pathway analyses. We bring together the modeling results and present synthetic spectra at the end of 3. Opportunities for future work and observations are addressed in Section 4 and we sum up the highlight of this study in Section 5.

The setup of 3D GCM simulations
We simulate the 3D global thermal and dynamical properties of HD 80606 b's atmosphere with the SPARC/MITgcm model (Showman et al. 2009;Adcroft et al. 2004). The SPARC/MITgcm model has been extensively applied to study the circulation of gas giants and brown dwarfs (e.g. Showman et al. 2009;Parmentier et al. 2016;Powell et al. 2019;Steinrueck et al. 2019;Tan & Komacek 2019). For radiative transfer, we assume chemical equilibrium and use the correlated-k method with 8 Gauss points within each of the 11 wavelength bins. A solar spectrum is used as an analog star for HD 80606. The simulations employ a resolution of C32 on the cubed-sphere grid and 53 vertical layers over 170 and 2.4×10 −6 bar. While the rotation period is one of the fundamental planetary properties that control the climate dynamics, it is challenging to place tight constraints on the rotation state of HD 80606 b. In this work, we follow Hut (1981Hut ( , 1982 and assume that the planet is "pseudo-synchronized" with its host star, i.e., the planet's rotation period is equal to the revolution period of a circular orbit at periastron, which is the dynamically stable configuration. The pseudo-synchronous rotation period is determined by the orbital period and eccentricity (Hut 1981) τ rot = τ orb × (1 + 3e 2 + 3 8 e 4 )(1 − e 2 ) 3/2 1 + 15 2 e 2 + 45 8 e 4 + 5 16 e 6 .
(1)  (2018) performed high-resolution spectral analysis and found the metallicity of HD 80606 to be about twice the solar value. Additionally, tidal dissipation is expected to increase the internal heating of eccentric planets (Leconte et al. 2010). Therefore, we further explore the effects of atmospheric metallicity and internal heating. We performed models for 1 time, 5 times solar metallicity, and internal temperature T int = 100 , 400 K, respectively. A summary of the model parameters used in this study are listed in Table  1. A periodic steady state is found after just a few orbits, where the temperature and wind structures at the same orbital position between different orbits show negligible changes. The output after nine orbits of simulation runs is used in this work.

The rise of temperature and winds around periastron
Owing to the high eccentricity of HD 80606 b, the stellar irradiation received by the planet shifts drastically from that similar to Earth at apoastron to about 800 times higher at periastron, raising the equilibrium temperature from ∼ 300 K to ∼ 1500 K (Figure 1). The planet spends more than half of the orbital period with a distance greater than 0.7 AU to its host star and only about two days with a distance less than 0.1 AU. We can in fact regard HD 80606 b as a warm Jupiter that enters a brief period of hot Jupiter regime every ∼ 100 days.
The thermal and dynamical responses to the rapid heating near periastron are illustrated in Figure 2. The pulse heating excites planetary-scale waves that associate with the momentum transport of the mean flow (Showman & Polvani 2011;Tsai et al. 2014;Showman et al. 2020). The patterns resemble those commonly seen in the spin-up stage for circulation under stationary forcing (Arnold et al. 2012;Debras et al. 2020;Hammond et al. 2020). Typically, the zonal jet on tidally locked planets takes about 10-100 Earth days to approach equilibrium from the initial rest state (i.e. the spin-up process; Showman & Polvani 2011;Hammond & Pierrehumbert 2018;Hammond et al. 2020). Hence, the rapid heating of ∼ 2 days near periastron does not provide enough time for the jet to fully develop before the forcing fades out. The jet speed evolves steadily throughout the orbit between 500 m/s and 2000 m/s but the circulation pattern outside of the periastron passage does not vary much and is similar to that shown at apoastron.
We can gain physical insight into the transition of dynamical regimes by comparing the circulation of HD 80606 b at different orbital phases with the grid of models with different stellar fluxes and rotational periods in Showman et al. (2015). The pseudosynchronous rotational period of HD 80606 b is close to the median rotation case (Ω med ; hereafter following the notation in Showman et al. (2015)) in Showman et al. (2015). The eccentric planet swings between the hot state (HΩ med ) and the cold state (CΩ med ). Showman et al. (2015) find that the hot state is dictated by the day-night thermal forcing while the cold state exhibits zonal symmetry and is mainly driven by rotation effects.
The high eccentricity makes HD 80606 b stay in HΩ med only for less than two days before returning to CΩ med . For most of the orbit, it shares circulation features similar to CΩ med , except the equatorial wind speed is faster than that in CΩ med . In this regime, the planet received weaker irradiation with the latitudinal temperature gradient predominating over the day-night temperature difference. The stellar irradiation swiftly increased by two orders of magnitude as the planet approaches periastron and enters HΩ med regime. In this regime, the thermal and wave structures resemble those in HΩ med but without the fast equatorial jet. This feature can be understood by considering the radiative timescale of an atmosphere with characteristic pressure P and equilibrium temperature T eq (e.g. Showman & Guillot 2002;Mayorga et al. 2021), where c P , g, and σ are heat capacity, gravity, and the Stefan-Boltzmann constant, respectively. τ rad of HD 80606 b is in the order of hours near periastron, which is consistent with the estimate from the phase curve observation (Laughlin et al. 2009;de Wit et al. 2016). Therefore, the thermal response of the atmosphere (∼hours) is much faster than the growth time of the jet (∼10-100 days). This thermal forcing initiates the Matsuno-Gill type (Matsuno 1966;Gill 1980;Showman & Polvani 2011) standing waves that reset the mean flow but at the same time pump momentum toward the equator. In the end, although the flash forcing does not allow the jet to develop during the short periastron passage, it contributes to accelerating the equatorial wind over time out of the periastron passage. This explains the stronger (∼ 1000 m/s faster) equatorial jet of HD 80606 b outside the periastron passage, while comparing to CΩ med in Showman et al. (2015) without the "flash forcing". The mean equatorial jet speed over a period is shown in Figure 3, where the sinusoidal cycle associates with wave-pumping acceleration for about half of the orbit and gradually slowing down of the jet for the second half. Figure 4 presents the overall response near periastron for different metallicities and internal temperatures, showing the dayside averaged temperature and root-mean-square (RMS) vertical wind as a function of time relative to periastron. A transient thermal inversion occurs immediately after periastron in all cases. Both temperature and vertical wind remain elevated until transit (∼ 5.6 days after periastron) with respect to those before periastron. Compared to the model with solar metallicity and T int = 100 K (referred as the nominal model hereafter), we find that a higher metallicity increases the atmospheric opacity and enhances the shortwave absorption, leading to a warmer stratosphere and a more pronounced updraft near periastron. On the other hand, the temperature and vertical wind in altitudes above the 1-bar level are not directly sensitive to internal heat.  90 species with 1028 thermochemical reactions and 60 photodissociation reactions. We adapted VULCAN to have temperature and stellar irradiation that continuously varied with time according to the orbital position of eccentric planets. The time series of GCM output provides the dayside average properties as input for VUL-CAN to simulate the compositional variation. HD 80606 is a G5 star with effective temperature of 5645 K and we used the solar spectrum (Gueymard 2018) as an analogue. Since the planet is assumed to be in pseudo-synchronous rotation, we define the dayside at periastron to be the dayside hemisphere throughout the orbit, i.e. our 1D photochemical model represents the dayside that is quasi-synchronous around periastron and keeps track of the same geographical hemisphere for the rest of the orbit. In this way, the dayside model captures the shifts of stellar irradiation near periastron, while the diurnal cycle of instellation makes the dayside and global average temperatures and wind almost identical away from periastron. We do not include the compositional feedback to the 3D GCM self-consistently and will discuss the limit of our 1D chemical model in 4.3. The same planetary parameters in Table 1 and a zenith angle of 58 • are used in the dayside-average photochemical model. Atmospheric mixing is commonly parameterized by eddy diffusion in 1D models. To track the variation of vertical mixing with time, we derive the eddy diffusion coefficients from the dayside averaged vertical wind (w rms ) at a given time according to the mixing length theory K zz = 0.1 w rms × H, with H being the pressure scale height and the factor of 0.1 following the empirical factor from tracer studies (e.g. Parmentier et al. 2013;Charnay et al. 2015). The temperatures and K zz are kept at the same values from those at the top boundary of the GCM (∼10 −5 bar ) when extending to 10 −8 bar for the photochemical model. Figure 5 shows the eddy diffusion profiles around an orbit, compared to the temperature scaling expression in Moses et al. (2021). Our wind-derived K zz shows a weaker pressure dependence, but the magnitude is consistent with the bounds given by the expression in Moses et al. (2021). We further explore the sensitivity to eddy diffusion in Section 3.3.3.
In this time-resolved photochemical model, one orbit of HD 80606 b is divided into 574 temporal grid points, with a separation of 1 hour near periastron to resolve the rapid variation and 13.3 hours for the rest of the orbit. The temperature and eddy-diffusion at each grid are updated from the GCM output and the UV flux is adjusted according to the orbital distance. The stepsize in the photochemical model is bounded by the grid size (e.g. 1 hour near periastron) to resolve the orbital evolution. To reduce the convergence time, we first run the chemical model at apoastron to the steady state. Then the steady-state composition at apoastron is used to initialize the timedependent model with orbital motion. With our initialization and the planet's long period (about 111 days), we find it takes less than 5 orbits to reach a periodic steady state. We have also tested the temperature dependence of UV cross sections of CO 2 and H 2 O ) and found negligible differences. After obtaining the composition, we use the open-source radiative-transfer tools HELIOS (Malik et al. 2019a) to generate emission spectra and PLATON (Zhang et al. 2019(Zhang et al. , 2020 for transmission spectra in Section 3.4.

Orbital-position independence of vertical quenching
First, following Visscher (2012), we will rely on equilibrium chemistry and timescale comparison to gain insight. Figure 6 illustrates the temperature profiles and the transitions of CH 4 -CO and NH 3 -N 2 , overlying the chemical timescales of CO and NH 3 . Thermochemical equilibrium predicts CH 4 dominates CO for the entire atmosphere at cold apoastron, while CO can take over CH 4 at altitudes above the 1-bar level at periastron as the temperature rapidly increased. Similarly, NH 3 is favoured over N 2 at apoastron between 1 bar and 0.1 mbar, while N 2 predominates at periastron. In reality, however, as the timescale of chemical conversion exceeds that of atmospheric transport at lower temperatures/pressures, the species retain uniform distributions above the levels where the two timescales are equal, i.e. the quench levels (e.g., Visscher et al. 2010;Tsai et al. 2018). The quenched abundances generally provide the basis distribution in the observable region (Moses 2014;Baxter et al. 2021;Fortney et al. 2020).
For eccentric planets, one might expect that the shift of temperatures and the strength of mixing would lead to temporal variations in the quenching process. Yet it is not the case for HD 80606 b. Figure  6 indicates that both CO and NH 3 are quenched in the deep atmosphere between 10 and 100 bar, where the radiative time scale is too long to manifest any temperature variations along the orbit. Hence, for a given K zz , the quench levels are independent of the thermal variations throughout the orbit. The elevated vertical wind ( Figure  2) after periastron also turns out to be not important for affecting the quench levels, because the steep decline of chemical timescales implies that the quench levels are not sensitive to merely a few factors of variation of K zz . Ultimately, we find the quenched abundances of the major carbon and nitrogen species remain unchanged throughout the orbit. In the next section, we will discuss how the chemical response is predominantly driven by photochemistry.

Photochemistry-driven response near periastron
Photochemistry produces abundant radicals, which can drive atmospheric composition to respond much faster (than what pure thermal kinetics allows) in the upper atmosphere. We now examine how the composition varies in response to the sudden shift of incident stellar flux and temperature during the periastron passage. The compositional response is depicted by the mixing ratio profiles as a function of time in Figure 7.
The stratosphere between ∼ 1 and 10 −4 bar quickly heats up when the planet is approaching periastron. It becomes warm enough (about half day before and one day after periastron) such that CO would have dominated over CH 4 above 1 bar in thermochemical equilibrium, as suggested by Figure 6 and confirmed in Figure 8. However, the duration of this hot state is orders of magnitude shorter than the CO-CH 4 chemical timescale ( Figure 6). We find photochemically produced radicals, such as H, OH, and S, to be key to promptly initiate the CH 4 -CO conversion. After periastron, the produced CO mixing ratio can exceed 10 −4 for a few days before getting dissipated by vertical mixing. An intriguing feature is that the slope of CO abundance in time in Figure 7 essentially captures the lifetime of the transient CO.
The variation of [CO]/[CH 4 ] ratio at 1 mbar is further illustrated in Figure 9.
[CO]/[CH 4 ] starts to rise just before the eclipse and peaks about 8 hours after periastron, where [CO]/[CH 4 ] 10 -10 4 , depending on the atmospheric metallicity and internal temperature. The JWST Cycle 1 observations probing the atmosphere above ∼ 10 mbar, especially by NIRSpec where CH 4 and CO have strong features within this wavelength range, will potentially provide the first real-time track of CH 4 -CO conversion.
We also find that the variation of H 2 O inversely follows CO, as oxygen is transferred from H 2 O to CO in the net reaction CH 4 + H 2 O − −− → CO + 3 H 2 . In addition, photochemical products, such as HCN and C 2 H 2 , are promoted near periastron and have significantly higher concentrations than their equilibrium abundances. The production of HCN follows the shared dissociation levels of CH 4 and NH 3 (Moses et al. 2011;Tsai et al. 2021) and dissipates similarly to CO after periastron. C 2 H 2 briefly reaches a high concentration of 10 −4 across the atmosphere (0.1 -10 −5 bar) around periastron. The transient peak of C 2 H 2 is a combination of C 2 H 2 being photochemically produced and thermodynamically favored at temperature 1000 K around periastron . Conversely, CO 2 momentarily falls back to the lower abundance in chemical equilibrium by about two orders of magnitude during this period. The peak of CO 2 in the upper atmosphere is restored just before transit, following the recovery H 2 O which produces OH radical turning CO into CO 2 . Lastly, while most substantial variation of the main carbon, oxygen, and nitrogen  Figure 6. The temperature profiles of our HD 80606 b model at apoastron (blue) and periastron (orange) overlaying the contours of the chemical timescales of CO and NH 3 , given by the expressions in Tsai et al. (2018). The dashed line shows where CH 4 and CO (NH 3 and N 2 ) are in equal abundance. The green circles indicate the quench levels corresponding to vertical mixing timescales of 10 4 -10 7 s taking into account mixing-length uncertainty (estimated from τ mix ∼ H 2 /K zz in the atmosphere below 1 bar). The quench levels are deep such that the temperature at the quench levels remains unchanged along the orbit. species occur near periastron, the overview of the variation of their abundances along the full orbit can be found in Figure A1.

The effect of orbit-induced quenching
To further understand how the eccentric orbit impacts the composition variations, we perform control experiments to isolate the effects of orbital motion. In this setup, we rerun the same photochemical model at each orbital position except keeping the planet fixed, allowing the composition to evolve until the steady state has reached. Figure 10 provides the mixing ratio profiles at the same orbital positions as Figure 7. Figure 10 shows that with enough time, CH 4 can fully convert to CO above ∼ 10 mbar within the synchronous-rotation period in the nominal case (about 1 day before and after periastron). The depletion of CH 4 subsequently restrains the production of HCN, making the variation of HCN completely differ from the continu- ous response in the nominal model with orbital motion. Outside of the synchronous-rotation period, CH 4 remains the dominant carbonbearing molecule. The elevated vertical mixing after periastron is somewhat more favorable for CH 4 than CO and subsequently suppresses CO 2 at steady state. This further confirms that in the nominal model with an eccentric orbit, the peak of CO 2 around transit is a result of the leftover CO.
We highlight the role of orbit-induced variations in Figure 11 by comparing the nominal model, fixed-orbit model, and chemical equilibrium for apoastron and periastron. While the difference between the equilibrium abundance and the fixed-orbit model informs 10 7 10 6 10 5 10 4 10 3 Mixing Ratio  us the effect of disequilibrium chemistry alone, the comparison between the nominal model and the fixed-orbit model allows us to tease out the orbital effects.
The close match between the nominal and fixed-orbit models in the left panel of Figure 11 indicates that the orbit-induced variations has minimal effect at apoastron, since the irradiation and temperature vary little away from periastron. In this cold state, CO maintained a quenched mixing ratio around 10 −5 , with a small amount of leftover from the hot state. It is predominantly vertical quenching, which is generally time independent (see the discussion in 3.2.1), that controls the major species with little contribution from the orbital motion. On the other hand, the atmosphere experiences rapid heating and a surge in the stellar flux over a few days during the periastron passage. The eccentric orbit motion takes over to dictate the thermal and chemical variations in this hot state. This can be seen in the right panel of Figure 11 that several main species undergo rapid change but do not have enough time reach steady state. CO production initiates around 1 mbar while the residual H 2 O and CH 4 from the cold state are partially destroyed. An important consequence is that under the intense UV flux, the excess CH 4 makes C 2 H 2 and HCN much more abundant above 10 mbar due to orbit-induced quenching.

The short-term and long-term response of sulfur species
We have discussed the behavior of carbon, oxygen, and nitrogen species so far. Next, we will examine the response of sulfur species. Figure 12 presents the dayside composition variation of several important sulfur species near periastron while that of the full orbit is shown in Figure 13. H 2 S is the thermodynamically stable sulfur bearing molecule in a H 2 atmosphere. Under shortwave irradiation where H atoms are produced in the upper atmosphere, H 2 S is attacked by H and branches to S or SH (Zahnle et al. 2016). The reaction that recycle sulfur back to H 2 S, has a high energy barrier and strong temperature dependence. For instance, the rate constant of Reaction (3) at about 0.1 bar increases from 10 −20 cm −3 s −1 to 10 −13 cm −3 s −1 as temperature climbs from 400 K to 1300 K near apoastron. Therefore, H 2 S can remain stable at higher altitudes with rising temperatures. The stable level of H 2 S tracks the thermal structure around periastron ( Figure 4) and remains elevated for a few days after periastron. The elevated abundance of H 2 S is important in bringing it to the H-rich photolysis region and initiating the production of elemental S , which sets up both short-term and long-term evolution along the orbit. For the short-term response during periastron passage, S is oxidized to SO and SO 2 in the upper atmosphere where OH is abundant, reaching a peak just before transit. Unlike other main species that respond primarily to the rapid irradiation near periastron and exhibit little variation for the rest of the orbit, sulfur also undergoes long-term cycles across the orbit, where it grows into elemental allotropes in the cold and reducing conditions. The evolution of sulfur species for a full orbit are illustrated in Figure 13. After periastron passage, the temperature gradually cools down and the lower UV flux allows the upper atmosphere to resume a more reducing condition. S starts to from allotropes in sequence: Once abundant atomic S is liberated from H 2 S by photochemistry, S 2 can form via either S + SH − −− → S 2 + H or S + H 2 S − −− → S 2 + H 2 . S 2 then continues with third-body assisted self reactions to form bigger allotropes (S x ). These growth steps following S − −− → S 2 , S 3 − −− → S 4 − −− → S 8 are captured in Figure 13. We find S 8 to be the main sulfur-bearing species in the gas phase above 0.1 bar outside the periastron passage. S 8 becomes saturated and condenses to form sulfur hazes between 1 and 0.01 bar with temperatures lower than about 360 K (Gao et al. 2017;Tsai et al. 2021). This long-term sulfur cycle is remotely analogous to the photochemical evolution after Shoemaker-Levy 9 impacts and shocks on Jupiter (Moses et al. 1995;Moses 1996).

The impact of sulfur on other species
Tsai et al. (2021) identify that sulfur can impact other species in a nonlinear way, including accelerating CH 4 -CO conversion with catalytic paths. Since CH 4 can only partially convert to CO near periastron, HD 80606 b provides an ideal setting to test the effects of sulfur on carbon conversion. To this end, we performed our nominal model except without including sulfur chemistry. Figure 14 compares the abundances from the models with and without sulfur six hours after periastron, when the CH 4 -CO conversion reaches maximum. Sulfur participation significantly enhances CO production between 10 −2 bar and 10 −4 bar, producing the peak of CO seen in Figure 7. We performed pathway analysis to look into the role of sulfur participation. Without sulfur, the main pathway for CH 4 -CO conversion is where the reaction involving water to form H 2 CO is the rate-limiting step. In the presence of sulfur, CH 4 goes through a catalytic pathway, where the rate-limiting step in (5), CH 2 + S − −− → HCS + H, is about 3 orders of magnitude faster than the reaction producing H 2 CO in (4) at 1 mbar. This reaction and pathway (5) rely on abundant S to be efficient, thus the peak of sulfur-enhanced conversion of CH 4 -CO follows S as shown in Figure 14.

Sensitivity to atmospheric metallicity and internal heat
We have seen that CH 4 , instead of CO, remains the dominant carbonbearing molecule for the majority of the orbit except for a few days during the periastron passage. Since CH 4 is less favored at high metallicity and high temperature, we now address whether the gen-eral behavior holds when we increase atmospheric metallicity and internal heat. We also test the sensitivity to vertical mixing.

Sensitivity to internal heat
The higher T int sets up a greater lapse rate in the deep layers below 10 bar, but the thermal and dynamical structure above 10 bar level remains insensitive to this increase in internal heat (Figure 4). We are aware that the GCM might take a much longer time to converge in the deep region (Carone et al. 2020;Mendonça 2020;Wang & Wordsworth 2020). Nevertheless, it is unlikely that the temperature and wind above 10 bar would be largely affected by the deep evolution. CH 4 remains the dominant carbon-bearing molecule despite the higher quenched CO abundance. Overall, our T int = 400 K models show qualitatively similar thermal and chemical variations to the same models but with T int = 100 K, as seen in Figure 9 (solid lines and dashed lines).

Sensitivity to metallicity
The temperature in the photosphere rises when metallicity increased (Figure 4), owing to the increased opacities (Drummond et al. 2018), mainly contributed by water. The combined thermal and chemical effects make the relative CH 4 and CO ratio more sensitive to the change of metallicity than raising the internal heat alone. Figure B1 compares the composition distributions with solar and five times solar metallicity at apoastron and periastron. When T int remained the same (upper panels), CO and CH 4 have very close quenched abundances when metallicity increased by 5 times. During the periastron passage, CO takes over CH 4 as the main carbon molecule above about 0.1 bar. When T int is raised to 400 K, the general trend remains the same but now the quenched CO abundance exceeds that of CH 4 . An intriguing feature is that in both of our 5 times solar metallicity models, a second CH 4 -CO conversion peak occurred deeper at about 0.1 bar due to stronger temperature inversion. As indicated in Figure 9, the maximum CO/CH 4 ratios with five times solar metallicity are 1-2 orders of magnitude higher than their solar-metallicity counterparts. Overall, the parameters we explored suggest that both supersolar metallicity and high T int are required to have CO as the predominant carbon molecule over CH 4 throughout the orbit. 3.0 × 10 18 Figure 13. The full-orbit evolution of the dayside composition of several sulfur species that form allotropes and hazes. S 8 particle mixing ratio is shown in n S8 /n gas , where n gas is the number density (cm −3 ) of the total gas. Note the elemental sulfur cycle consisted of S-S 2 -S 4 -S 8 sequence that spanned the orbit.

Sensitivity to K zz
Regarding the parameterizing uncertainty of vertical mixing, we performed sensitivity tests with varying eddy diffusion. Figure B2 compares the models with the nominal eddy diffusion coefficient derived from the GCM's vertical wind to those scaled by 3 and 1/3. We find the quench levels remain not too sensitive to shifting K zz . However, eddy diffusion can affect composition distributions in two major ways. First, vertical mixing controls the level where molecular diffusion becomes dominated and the species started to stratify following their own scale heights, i.e. homopause (Leovy 1982;Sinclair et al. 2020). This stratification effect is seen in the decay of the water profiles in Figure B2. Second, eddy diffusion controls the photochemical lifetime of the species in the stratosphere produced during the periastron passage. Stronger vertical mixing transports and dissipates these species more efficiently, resulting in the K zz dependence of the surplus CO at apoastron in Figure B2.

Synthetic Spectra
In this section, we present the synthetic spectra generated by our 1D dayside average model. We do not take into account the change of Earth-facing hemisphere and the contribution from the nightside but focus on the shape and features of the spectra. Figure 15 compares the emission spectra of HD 80606 b six hours before, eight hours after, and at the secondary eclipse. We find that the emission flux is sensitive to the temperature increase with higher metallicity a few hours before and at the eclipse. During this period, the eclipse and pre-eclipse emission spectra are ideal for constraining the atmospheric metallicity. After the eclipse, the temperatures from models with different metallicities start to converge as the strong temperature inversion dominates the flux, leading to similar emission spectra. However, while the Spitzer measurement at 4.5 micron favors a slightly enhanced metallicity between solar and 5 × solar values, none of our models are able to fit the Spitzer IRAC Ch1 and Ch2 data simultaneously. Figure 16 illustrates how the emission spectra of our nominal model evolved before and after the eclipse. Since the chemical conversion does not occur instantaneously, significant differences in the CH 4 and H 2 O band exhibit when assuming chemical equilibrium. We find that the spectra before the eclipse are generally characterized by prominent emission features and those after the eclipse start to show absorption features. At the eclipse, the smaller temperature gradient makes the spectra closer to blackbody emission. The transition can be easily seen by the strong features of CH 4 at 3-4 and 7-9 micron. Similar transitions are also found for H 2 O around 11 micron and HCN around 14 micron. Only minor features of CO at 4-5 micron, owing to the pressure level where CO absorbs the most coincide with the isothermal part of the atmosphere. Nevertheless, the CH 4 -CO conversion can potentially be inferred by tracking CH 4 and H 2 O (also see Figure C1 for spectral variation due to trace composition isolated from the thermal variation). Within this JWST cycle 1 window between 10 hours before and after the eclipse, CO 2 falls back to the lower equilibrium composition and cannot be detected in the emission spectra.
Lastly, Figure 17 shows the theoretical transmission spectra of HD 80606 b, where the transit occurs about 5.6 days after periastron. CO, CO 2 , and HCN exhibit significant disequilibrium features at 4-5 and 14 micron, respectively. The effects of orbit-induced quenching can also be seen by comparing the nominal kinetics model with the "fixed-orbit" model. Notably, leftover CO is carried over and contributes to the production of CO 2 and HCN, where the feature of HCN at 14 micron shows the most pronounced differences. Since CO 2 reached a peak value in the upper atmosphere during the transit, we suggest future transit observations to look for CO 2 to put our photochemical modeling to the test. Although SO and SO 2 also reached highest abundances before the transit, SO 2 does not reach high enough abundances to be detected in all models we explored 2 .

Uncertainties in sulfur kinetics
Considerable uncertainties exist in sulfur kinetics, especially the reactions involving reduced sulfur species (Zahnle et al. 2016;Tsai et al. 2021). One important initial step in forming sulfur allotropes is to make S 2 from SH or S. For the three-body recombination reaction that form S 2 from S, S + S M − −− → S 2 , we followed previous works and adopted the latest calculation from Du et al. (2008). This rate coefficient by Du et al. (2008) is consistent with Fair & Thrush (1969) but about three orders of magnitude smaller than another experimental rate Nicholas et al. (1979). Thus, we consider our S x formation in this study as a conservative limit. Another unsettled reaction that might be more relevant for HD 80606 b is S + CO M − −− → OCS, whose reaction rates at low temperatures spread across a few orders of magnitude in the literature. Tsai et al. (2021) found that this reaction is key to OCS formation and S 2 production on directly imaged gas giants, while Ranjan et al. (2020) have shown that it can influence the CO abundance in a N 2 -CO 2 atmosphere as well. We reiterate the need of further investigation to resolve the OCS recombination rate.

Photochemical Hazes
The combination of CH 4 carried over from the cold phase and the intense UV radiation seemingly entails an ideal condition for photochemical hazes near periastron. Indeed, our nominal models find copious C 2 H 2 (∼ 100 ppm) along with some hydrocarbon haze precursors, such as C 6 H 6 (peaked around 0.01-0.1 ppm). If monomer formation and growth occur fast enough (Seinfeld & Pandis 2016;Ohno & Okuzumi 2017), organic hazes could affect the periastron passage. We will leave simulating the organic aerosols variations for future work (Ohno et al. in prep). On the other hand, S 8 is produced outside of the periastron passage and condenses between 0.01 and 1 bar once the temperature cools down. S 8 hazes scatter strongly and reduce the emission in the infrared. Figure 18 illustrates 1 µm-sized S 8 hazes substantially lower the emission flux. Our models suggest that organic hazes might be promoted during the periastron passage and deeper sulfur hazes could prevail a large portion of the orbit.
Additionally, Na 2 S sodium sulfide could form via 2 Na + H 2 S − −− → Na 2 S(s) + H 2 , where Na 2 S particles condense out between 10 and 0.1 bar Mayorga et al. 2021). The 2 The opacity of the more abundant SO is not included due to lack of data. Na 2 S clouds can potentially suppress the emission flux but are evaporated and dissipated near periastron. In respect of the sulfur inventory, since sulfur is about 10 times more abundant than sodium in solar composition, it is unlikely that Na 2 S would substantially deplete sulfur.

Future modeling and observational implications
We applied a 1D photochemical model based on averaged 3D GCM output to simulate the time-varying chemistry on HD 80606 b for this study. In reality, there will be radiative feedback from the compositional variation (e.g., Drummond et al. 2016;Lavvas & Arfaux 2021) and 3D circulation can complicate the global distribution of chemical species (e.g., Mendonça et al. 2018;Drummond et al. 2020;Steinrueck et al. 2021). To assess the limit of our 1D photochemical model, we have compared the dynamical timescales and found that the horizontal transport timescale is generally shorter or comparable to the vertical mixing timescale. Since vertical quenching already occurred at 10 -100 bar, we expect horizontal transport simply lead to a uniform composition based on the vertically quenched abundance on the dayside. However, as shown in the study, photochemistry proceeds faster in the upper atmosphere and can have compelling interaction with the circulation. We encourage future work to treat the chemical feedback on the thermal structure consistently and/or utilize multidimensional photochemistry models (e.g., Venot et al. 2020;Tsai et al. 2021b;Baeyens et al. 2022) to explore the global properties of the time-varying, eccentric systems.
We have demonstrated that the quenched abundances of CO and CH 4 before the CH 4 -CO conversion kicks off around periastron are independent of the orbital phase. The measurement at the start of JWST cycle 1 observation should provide a baseline of the atmospheric properties. Our results suggest that [CO]/[CH 4 ] ratio is highly sensitive to metallicity and mildly sensitive to internal heat. Future spectral observations outside of the periastron passage of HD 80606 b can help us to place tighter constraints. Since CO 2 is typically more sensitive to atmospheric metallicity, future transit observations will improve our ability to constrain atmospheric metallicity and validate the photochemical process.
We found that sulfur species are oxidized after periastron and SO and SO 2 reach peak abundances near transit, while sulfur allotropes (S x ) exhibit a long-term growing sequence along the orbit. SH and S 2 are strong absorbers in the optical and NUV (e.g., S 2 emission has been reported after the impact of Shoemaker-Levy 9 on Jupiter (Noll et al. 1995;Zahnle et al. 2009)). S 3 absorbs effectively near 400 nm and S 4 has a weaker, broad-band absorption between 450 and 600 nm, whereas S 8 has stretching mode features around 20 micron. We propose future observations sampling a wider coverage of the full phase to track the sequence of polysulfur species.
Finally, although we focused on HD 80606 b in this study, it would be illuminating to learn whether our findings can be generalized to other eccentric systems. While detailed atmosphere variation depends on the interplay between radiative and chemical processes, we expect photochemistry to play a main role in the time-dependent chemistry on other eccentric systems, given its short timescale. Exploring the key parameter space such as eccentricity and the irradiation temperature at periastron with coupling radiative transfer model and photochemical kinetics would be an ideal setting to address this question. 6 hrs to eclipse 6 hrs to eclipse in chemical eq. +8 hrs to eclipse +8 hrs to eclipse in chemical eq.

CONCLUSIONS
We present an integrated modeling framework to study the atmospheric responses to the irradiation variation of the highly eccentric gas giant HD 80606 b. We applied a 3D GCM to simulate the fast-evolving climate dynamics and fed the dayside average atmosphere profiles to a time-dependent 1D photochemical model for the response of chemical compositions. We made predictions with synthetic spectra and demonstrate that our knowledge of atmospheric physics and chemistry will soon be put to the test with JWST observation.
We summarize our findings as follows, • The climate of HD 80606 b oscillates between the cold and hot states that are governed by distinct and well-studied dominating processes (day-night forcing vs. baroclinic instability).
• Vertical quenching occurs in the deep atmosphere and is generally independent of the orbital position. Photochemistry from the intense UV radiation near periastron is the main driver of the compositional variation.
• CH 4 is partly converted to CO during the periastron passage and can be tracked by the JWST Cycle 1 GO observation. We have isolated the orbital effect in the model and showed that orbit-induced quenching plays an important role in this eccentric system.
• We found CO 2 reached a peak abundance near transit and could potentially be detected in transmission. Sulfur exhibits short-term response with oxidized species and long-term variations with sequence of allotropes cycle, where S 8 is the main sulfur-bearing molecule above 0.1 bar away from periastron.
• Our photochemical results suggest two classes of hazes as S 8 condense out for most of the orbit while organic haze precursors are promoted near periastron. In this appendix, we show the variation of the abundances of several important carbon, oxygen, and nitrogen species along the full orbit to reiterate that the main chemical response occur during the short periastron passage while the composition of these species remain invariant for most of the orbit. Here, we present the composition distributions at apoastron and periastron when we vary the atmospheric metallicity, internal heat, and vertical mixing.  H 2 O HCN 6 hrs to eclipse 6 hrs to eclipse with swapped compostion +8 hrs to eclipse +8 hrs to eclipse with swapped compostion Figure C1. Emission spectra of HD 80606 b six hours before the eclipse and eight hours after the eclipse, similar to Figure 16, but here we show the spectra with swapped composition in grey (temperature profile -6 hr before the eclipse with composition +8hr after the eclipse and vice versa) to demonstrate the spectral changes in chemical composition isolated from thermal changes.