The effect of internal gravity waves on cloud evolution in sub-stellar atmospheres
Amy Parent, Ruth E. Falconer, Elspeth K. H. Lee, Karen A. Meyer, Craig R. Stark
AAstronomy & Astrophysics manuscript no. parent_2019 c (cid:13)
ESO 2020February 28, 2020
The effect of internal gravity waves on cloud evolution insub-stellar atmospheres
A. Parent , R. E. Falconer , G. K. H. Lee , K. A. Meyer , and C. R. Stark Division of Computing and Mathematics, Abertay University, Kydd Building, Dundee DD1 1HG, UKe-mail: Mathematics, School of Science & Engineering, University of Dundee, Nethergate, Dundee DD1 4HN, U.K. Atmospheric, Oceanic and Planetary Physics, Department of Physics, University of Oxford, Oxford OX1 3PU, UKReceived ... / Accepted ...
ABSTRACT
Context.
Sub-stellar objects exhibit photometric variability, which is believed to be caused by a number of processes, such asmagnetically-driven spots or inhomogeneous cloud coverage. Recent sub-stellar models have shown that turbulent flows and waves,including internal gravity waves, may play an important role in cloud evolution.
Aims.
The aim of this paper is to investigate the e ff ect of internal gravity waves on dust nucleation and dust growth, and whetherobservations of the resulting cloud structures could be used to recover atmospheric density information. Methods.
For a simplified atmosphere in two dimensions, we numerically solved the governing fluid equations to simulate the e ff ecton dust nucleation and mantle growth as a result of the passage of an internal gravity wave. Furthermore, we derived an expressionthat relates the properties of the wave-induced cloud structures to observable parameters in order to deduce the atmospheric density. Results.
Numerical simulations show that the density, pressure, and temperature variations caused by gravity waves lead to an increaseof the dust nucleation rate by up to a factor 20, and an increase of the dust mantle growth rate by up to a factor 1 .
6, comparedto their equilibrium values. Through an exploration of the wider sub-stellar parameter space, we show that in absolute terms, theincrease in dust nucleation due to internal gravity waves is stronger in cooler (T dwarfs) and TiO -rich sub-stellar atmospheres. Therelative increase, however, is greater in warm (L dwarf) and TiO -poor atmospheres due to conditions that are less suited for e ffi cientnucleation at equilibrium. These variations lead to banded areas in which dust formation is much more pronounced, similar to thecloud structures observed on Earth. Conclusions.
We show that internal gravity waves propagating in the atmosphere of sub-stellar objects can produce banded cloudsstructures similar to that observed on Earth. We propose a method with which potential observations of banded clouds could be usedto estimate the atmospheric density of sub-stellar objects.
Key words. brown dwarfs – stars: atmospheres – hydrodynamics – waves
1. Introduction
Brown dwarfs are low-mass, sub-stellar objects below the hydro-gen burning limit, with masses between 13 M Jup and 70 M Jup . Asa consequence, their atmospheres are su ffi ciently cool for the for-mation of dust clouds. The clouds are observed through their ef-fect on spectral features and their association with infrared spec-troscopic variability, which is believed to be caused by patchyclouds. Numerous observations (see for example Buenzli et al.2014) show that a large portion of known brown dwarfs exhibitphotometric variability. According to Biller (2017), over 10% ofknown brown dwarfs show a variation of 1% or more, and over50% exhibit a variation of 0.1% to 0.5% or more.Explaining spectral variability is necessary to understandingthe L / T transition (Vos et al. 2019). The L and T componentsof the Luhman 16AB system, for example, show vastly di ff er-ent patterns (Gillon et al. 2013): Luhman 16B, a T dwarf, ex-hibits strong, fast-changing periodic spectral variations, whilethe L dwarf Luhman 16A exhibits no periodic pattern. A modelfrom Saumon & Marley (2008) and Marley et al. (2010) pro-poses the sinking of parts of the cloud deck, creating thinner,patchy clouds as an explanation for the change in variability pat-terns around the L / T transition. A study by Stark et al. (2015)proposes the electrostatic disruption of charged cloud particles as a mechanism through which inhomogeneous coverage couldbe caused. While inhomogeneous dust cloud coverage (Helling& Casewell 2014) is believed to be the main cause for browndwarf variability, other theories, such as temperature variations(Robinson & Marley 2014) or fingering convection (Tremblinet al. 2016), propose cloud-free models as an explanation.Internal gravity waves have been simulated in main sequencestars (Alvan, L. et al. 2014). Internal gravity waves triggeredby fingering convection have been modelled in objects rang-ing from main-sequence stars to brown dwarfs (Garaud et al.2015), and they can reach wavelengths much larger than thesource perturbation. Simulations of atmosphere patches by Frey-tag et al. (2010) show that internal gravity waves are also presentin sub-stellar atmospheres, triggered by downdrafts caused byconvection patterns, and they are theorised to be one of the mainphenomena responsible for transporting dust in the upper atmo-spheric layers.Further characterising the impact of gravity waves on dustcloud evolution can advance the understanding of cloud struc-tures in brown dwarfs. Helling et al. (2001) showed that higherfrequency acoustic waves, triggered by turbulent flow in browndwarfs atmospheres, can have a strong impact on cloud forma-tion: by carrying lower temperature perturbations, the passage of
Article number, page 1 of 10 a r X i v : . [ a s t r o - ph . S R ] F e b & A proofs: manuscript no. parent_2019 waves can temporarily create favourable conditions for dust nu-cleation in otherwise dust-hostile environments, leading to theformation of dust over time.Internal gravity waves are a type of fluid wave that occursin atmospheres and oceans. Their defining characteristic is thatgravity, in the form of buoyancy, is the restoring force that al-lows disturbances to propagate. Internal gravity waves can be ob-served in the Earth’s atmosphere through their e ff ect on clouds.On Earth, a wave cloud is formed from the passage of an internalgravity wave, triggered by stable air flowing over relief. The ver-tical displacement of the air forces it to oscillate as the buoyancyforce tries to restore equilibrium. As the wave propagates, at thewave peaks the displaced air rises and cools resulting in watervapour condensing, forming droplets and clouds; at the wavetroughs, the clouds evaporate due to adiabatic heating, leadingto clouds that have a distinct banded structure.In the case of a gas giant planet or brown dwarf, dust cloudsare formed instead of water clouds but an analogous process canoccur. In this context, wave clouds can be induced as a result ofexternal fluid motion triggering turbulent flow (such as fingeringconvection in deeper layers of the atmosphere), whereas gas flowover relief would be expected to be the main cause in the case ofa rocky terrestrial exoplanet (Roeten et al. 2019). In a sub-stellaratmosphere, oscillating parcels of gas can trigger the nucleationof seed particles and enhanced surface mantle growth, formingbanded cloud structures. The nucleation rate is a function of den-sity of the nucleating species and the atmospheric temperature.If the passage of the internal gravity wave perturbs the local ther-modynamic structure of the atmosphere it can give enhanced nu-cleation in localised regions.The aim of this paper is to investigate and characterise thee ff ect of internal gravity waves on the evolution of dust cloudsin the atmospheres of sub-stellar objects and its consequencesfor cloud variability. This paper presents a novel mechanism forpotentially diagnosing the gas density of sub-stellar atmospheresfrom observations of the resulting cloud structures formed fromthe passage of internal gravity waves. In Sect. 2 the basic atmo-spheric model of internal gravity waves, nucleation and mantlegrowth is described; in Sect 3 the numerical methods used tosimulate internal gravity waves are presented; in Sect. 4 the re-sults of the simulations are presented and discussed; Sect. 5 sum-marises and discusses the consequences of the results includinga possible way of connecting observations to the wave dispersionrelation to diagnose the atmospheric density.
2. Sub-stellar internal gravity waves
For a vertical slice of a sub-stellar atmosphere in hydrostaticequilibrium, the coupled equations of fluid dynamics governingthe evolution of the fluid velocity u , the fluid density ρ , and thepressure p , of an atmospheric parcel, under the influence of grav-ity g are: ∂ρ∂ t + ∇ · ( ρ u ) = , (1) ρ (cid:34) ∂ u ∂ t + ( u · ∇ ) u (cid:35) = −∇ p − ρ g , (2) p ρ − γ a = const. , (3)where γ a is the ratio of specific heats (for a diatomic gas γ = / u = ∂/∂ t =
0, giving the fol-lowing equilibrium relationship: ∇ p = − ρ g , (4) where the subscript ‘0’ denotes an equilibrium quantity. For sim-plicity, in order to capture the fundamental physics, this paperfocuses on the e ff ect of gravity waves in the linear regime. Wecan linearise Eqs. (1)-(3) by decomposing each variable Q intoits equilibrium and perturbed value, so that Q = Q + Q . Inthe non-linear regime Q (cid:29) Q , and powers of Q higher than1 can be discarded. For clarity, the subscript ‘1’ is omitted infurther equations. Linearisation yields the final system of fluidequations: ∂ρ∂ t = −∇ · ( ρ u ) , (5) ∂ u ∂ t = − ∇ p ρ − ρ g ρ . (6)To model internal gravity waves, where we deal with incom-pressible flows ∇ · u =
0, we adopted a vorticity-stream functionformulation by introducing the vorticity ζ and stream function ψ defined by, ζ = ∇ × u , (7) u = ∇ × ψ . (8)Therefore, the governing fluid equations for incompressibleflows in the linear regime become, ∂ρ∂ t = − ( ∇ × ψ ) · ∇ ρ , (9) ∂ ζ ∂ t = −∇ × (cid:32) ρ g ρ (cid:33) , (10) ζ = ∇ × ( ∇ × ψ ) . (11)where the baroclinic term ( ∇ ρ × ∇ p /ρ ) vanishes since the prop-agation of internal gravity waves is considered to be an adiabaticprocess. In an atmospheric vertical plane ( x , y ): u = (cid:16) u x , u y , (cid:17) = (cid:32) ∂ψ∂y , − ∂ψ∂ x , (cid:33) , (12) ψ = (0 , , ψ ) , (13) ζ = (cid:32) , , ∂ ψ∂ x + ∂ ψ∂y (cid:33) , (14) g = (0 , g, . (15)This yields the system of equations ∂ζ∂ t = − gρ ∂ρ∂ x , (16) ∂ρ∂ t = ∂ρ ∂y ∂ψ∂ x , (17) ζ = −∇ ψ, (18)where ∇ = (cid:32) ∂ ∂ x + ∂ ∂y (cid:33) . (19)To obtain the gravity waves’ dispersion relation, Eq. (18) isderived with respect to time, and Eq. (16) substituted for ∂ζ/∂ t : ∂ζ∂ t = − ∇ (cid:34) ∂ψ∂ t (cid:35) = − gρ ∂ρ∂ x . (20) Article number, page 2 of 10. Parent et al.: Sub-stellar internal gravity waves Di ff erentiating Eq. (20) with respect to time, and then sub-stituting Eq. (17) for ∂ρ/∂ t yields: ∇ (cid:34) ∂ ψ∂ t (cid:35) = gρ ∂∂ x (cid:32) ∂ψ∂ x ∂ρ ∂y (cid:33) . (21)At equilibrium, ρ does not vary along x , Therefore: ∇ (cid:34) ∂ ψ∂ t (cid:35) = gρ ∂ρ ∂y ∂ ψ∂ x = − N ∂ ψ∂ x , (22)where N = (cid:115) − gρ ∂ρ ∂y , (23)is the Brunt-Väisälä buoyancy frequency. In the case of a uni-formly stratified atmosphere, assuming a solution of the form ψ ≈ exp [ − i ( ω t + k x x + k y y )] the dispersion relation for internalgravity waves becomes (Sutherland 2010; Vallis 2017), ω = (cid:115) N k x k x + k y , (24) = N cos θ. (25)Therefore, when the atmospheric density (or equivalently ve-locity) is perturbed, corresponding oscillations in ρ , ψ , and ζ aretriggered, that occur at the Brunt-Väisälä buoyancy frequency N . As the wave propagates through the atmosphere, the densityvariations can a ff ect the resulting nucleation and mantle growthrates. To quantify the impact of passing waves on dust formation, weused the equation of modified classical nucleation theory pre-sented by Gail et al. (1984); Helling et al. (2001), which definesthe nucleation rate J ∗ , the number of nucleating centres formedper second per unit volume [m − s − ], as: J ∗ = n x τ Z exp (cid:34) ( N ∗ −
1) ln S − (cid:18) T θ T (cid:19) N ∗ − N ∗ − / (cid:35) , (26)where T is the temperature; τ is the seed growth time scale forthe gaseous nucleation species x ; N ∗ is the size of the criticalcluster; n x is the number density of the nucleating species; Z isthe Zeldovich factor; S is the supersaturation ratio; defined asfollows: τ = n x v rel , x N / ∗ A , (27) N ∗ = + (cid:32) T θ T ln S (cid:33) , (28) Z = (cid:20) T θ π T ( N ∗ − / (cid:21) / , (29) S = p x p sat , x , (30) where,ln ( p sat , x ) = . − . T [dyn] , T ∈ [500 , , (31) v rel , x ≈ (cid:114) k B T π m x , (32) T θ = π r σ k B , (33) r = (cid:32) Am p πρ m (cid:33) / , (34)and p x is the partial pressure of the nucleating species x ; p sat , x isthe saturation vapour pressure of the nucleating species x ; σ isthe surface tension of the nucleating species; r is the hypothet-ical monomer radius; A = π r is the hypothetical monomersurface area; m x is the mass of a monomer particle; m p is theproton mass; A is the atomic weight of the monomer; ρ m is thedust material density; and v rel is the thermodynamic equilibriumvelocity for the nucleating species studied (TiO for this paper).For this paper, we assumed a constant value for the surfacetension of TiO , σ TiO = .
618 J m − (see Helling et al. (2001);Lee et al. (2015)), and a density ρ m ≈ − for TiO .In the context of an internal gravity wave propagating through asub-stellar atmosphere, the wave perturbs the local atmosphericgas density and hence the temperature (via T ρ − γ a = const.) in anadiabatic process. As a result, the passage of the wave perturbsthe nucleation rate J ∗ . Once nucleation has established a material surface onto whichmaterial can accumulate, dust growth occurs via gas-phase sur-face chemistry (Eq. (24) in Helling & Woitke (2006)). We con-sider a spherical dust grain of radius a , of mass m d , and let n bethe number density of the gas phase. The dust grain absorbs gasmolecules at a rate απ a n (cid:104) v (cid:105) , where (cid:104) v (cid:105) is the mean gas molec-ular speed and α is the sticking probability that a molecule isabsorbed (the sticking factor). Therefore, the mass of the dustgrain m d evolves in time asd m d d t = απ a nm (cid:104) v (cid:105) = απ a ρ (cid:104) v (cid:105) , (35)where ρ = nm is the gas mass density. The mass of a dust graincan be written as m d = π a ρ m , where ρ m is assumed to be con-stant. Therefore, the time evolution of the radius of a dust grainisd a d t = αρ (cid:104) v (cid:105) ρ m = γ, (36)where γ is the growth rate from absorption in units of [ms − ].Eq. (36) is the archetypal equation describing the absorption ofmaterial onto the surface of a dust grain. It is consistent with thedust growth equations presented in Helling et al. (2001), albeit ina much simplified form but still encapsulating the fundamentalunderlying physics. Furthermore, Eq. (36) is also consistent withmantle growth via ion accretion when dust grains are immersedin a plasma (Eq. (18) in Stark & Diver (2018)). Without loss ofgenerality, to investigate the e ff ect of internal gravity waves onthe mantle growth rate we simplify the expression by introducing ρ s , the density of the gas-phase accreting species, defined as ρ s = Article number, page 3 of 10 & A proofs: manuscript no. parent_2019 f s ρ , where 0 ≤ f s ≤ γ = f s ρv rel , s ρ m , (37)where we have set α = (cid:104) v (cid:105) = v rel , s . Introducing f s allows us to generalise the e ff ect ofdi ff erent gas-phase species, with varying relative abundances inthe gas-phase, participating in surface chemistry leading to man-tle growth.
3. Numerical simulations
The linearised governing equations can be cast in non-dimensional form, defining τ = t / T , (38) ξ = x / L , (39) λ = y/ L , (40) a = ζ T , (41) b = ρ L / M , (42) c = ψ T / L , (43) β = g T / L , (44)where L , T , and M are characteristic values of length, time,and mass respectively. Therefore, in non-dimensional form,Eqs. (16)-(18) become ∂ a ∂τ = − β a ∂ b ∂ξ , (45) ∂ b ∂τ = ∂ a o ∂λ ∂ c ∂ξ , (46) a = − ∂ c ∂ξ − ∂ c ∂λ . (47)Casting the model equations in non-dimensional form letsus observe the characteristic behaviour of internal gravity waveswithout loss of generality. To solve the system of fluid equations (45) - (47) numerically,we used a combination of the leapfrog method for Eq. (45) andEq. (46), and Successive Over-Relaxation (SOR) for Eq. (47)(Vetterling et al. 1992; Mittal 2014). SOR requires that the valueof ζ from Eq. (18) be known at the boundaries of the domain. Inorder to minimise the artefacts caused by boundaries, we used adomain large enough that over short timescales, the waves’ per-turbations do not reach the boundaries. Additionally, values of ρ , ψ , ζ , and their spatial first-order derivatives were interpolatedusing third-order polynomial interpolation. The internal gravitywaves were triggered by creating a Gaussian density perturba-tion initial condition in the centre of the numerical domain, thatwas modulated in time by a sine wave with a period equal to amultiple of the local buoyancy frequency, ρ = ρ A sin ( ωτ ) exp (cid:16) − ς r (cid:17) , (48)where ρ A is the maximum amplitude of the perturbation; ς is thespread parameter; and r is the distance to the centre of the numer-ical domain. The resulting wave solutions were used to calculate Table 1.
Parameters used for numerical simulations
Parameter Value n x xn y y h × − spacial mesh increment ω [0 . , .
0] driving perturbation frequency n τ
50 total no. of time steps d τ T ω / n τ temporal mesh increment Table 2.
Models presented in Figs. 1 and 3
Key T e ff log g descriptionBD 1500 K 5 . . . . . The aim of this paper was to investigate the e ff ects of internalgravity waves on the evolution of dust clouds in sub-stellar atmo-spheres. To this end, we considered a brown dwarf atmospherecharacterised by T e ff = g = . rift -P hoenix model atmosphere and cloud forma-tion code (Hauschildt & Baron 1999; Helling et al. 2004; Helling& Woitke 2006; Witte et al. 2009, 2011). We note that the atmo-spheric extension y (Fig. 1; middle panel) is measured from thetop of the atmosphere.The atmosphere model used as input is one-dimensional;for this study, we assumed a horizontally uniform atmosphereat equilibrium, and expanded the model to two dimensions. Wechose the characteristic length L (see Sect. 3.2) to be the heightof the simulation domain (10 m). The numerical domain wasdefined by the geometric extension y ∈ [25 km ,
35 km]; the pres-sure p ∈ [3 × − bar , × − bar]; the temperature T ∈ [700 K ,
850 K]; and the density ρ ∈ [10 − kg m − , − kg m − ].For these conditions, n TiO / n = const. ≈ − (see Fig. 2). Theperturbation length scale was of the order of 10 m.
4. Results
In a stratified sub-stellar atmosphere, perturbing the backgrounddensity by vertically displacing a fluid parcel triggers vertical os-cillations as the buoyancy force tries to restore equilibrium. Theresulting density variations propagate through the atmosphereat the Brunt-Väisälä buoyancy frequency. Figure 3 shows thebuoyancy period T N = π/ N as a function of pressure p forthe sub-stellar atmosphere characterised by T e ff = g = .
0. The buoyancy period is of the order of 10 s acrossthe extent of the atmosphere, with the maximum period occur-ring at high atmospheric pressures. Also shown in Fig. 3, forcontext, are the buoyancy periods profiles of typical L and T
Article number, page 4 of 10. Parent et al.: Sub-stellar internal gravity waves T / K BDTDLDLGTDLGLD10 y / k m p /bar / kg m Fig. 1.
Atmospheric diagram ( p , T ; top panel), ( p , y ; middle panel),and ( p , ρ ; bottom panel) for the D rift -P hoenix brown dwarf atmospheremodel used in numerical simulations (BD: T e ff = , log g = . g ) = .
5) and low-gravity(LGLD, LGTD: log( g ) = .
0) L dwarf ( T e ff = T e ff = dwarfs (periods of the order of 10 s to 100 s), and low-gravity Land T dwarfs (longer buoyancy periods, of the order of 1000 s).In comparison, Buenzli et al. (2014) observe spectroscopic vari-ations on timescales of 100 s to 1000 s. The key parameter distin-guishing between the buoyancy period for di ff erent atmosphericmodels is the surface gravity; since T N ∝ g − / , objects withlower surface gravity will have a longer buoyancy periods. In thecase of a neutrally stratified atmosphere (i.e. N = p /bar n T i O / n Fig. 2.
Plot of the partial number density of TiO at equilibrium for themodel used in simulations (log g = . , T e ff = − bar), n TiO / n is constant around10 − . frequency of the wave is set by the frequency of the source forfrequencies below the buoyancy frequency. In sub-stellar atmo-spheres the length scale of convective motions deep in the atmo-sphere is related to the atmospheric pressure scale height by afactor between 1 and 10 (Marley & Robinson 2015; Tremblinet al. 2019) – for example, l conv ≈ km in Freytag et al. (2010)– varying with timescales of the order of 10 s (Tremblin et al.2019). p /bar T N / s BDTDLDLGTDLGLD
Fig. 3.
Buoyancy period T N = π/ N for the D rift -P hoenix brown dwarf(log g = . , T e ff = s to10 s) is consistent with the atmospheric models and simulation resultspublished by Freytag et al. (2010). In lower-density brown dwarfs, thebuoyancy period rises to the order of 10 s. Figure 4 shows the characteristic St Andrews cross patternemanating from density perturbation located at the centre of thenumerical domain, for the example driving frequency ω = . N and a driving perturbation amplitude ρ A = . ρ . In an adiabaticprocess the density variations are accompanied by sympatheticvariations in temperature ( T ∝ ρ γ a − = ρ / ) that go on to a ff ectthe local nucleation rate and surface mantle growth rate. Figure 4shows that areas where ρ is reduced by the passage of internalgravity waves are the areas with the largest increase of dust nu-cleation. The number density n x of the nucleating species is a Article number, page 5 of 10 & A proofs: manuscript no. parent_2019 component of the total gas density ρ = (cid:80) s m s n s (See Fig. 2,showing the relationship between n and n TiO ). Variations of ρ propagated by passing gravity waves are therefore reflected inthe partial density of the nucleating species n TiO , which in turnsleads to an adiabatic change in temperature both of the total gasphase and its components. The primary e ff ect that drives thevariation in nucleation rate is the temperature variation: as thedensity increases (decreases) there is a corresponding adiabaticincrease (decrease) in the temperature. Increasing (decreasing)the temperature of the nucleating species in this way, decreases(increases) the supersaturation ratio, S . E ffi cient nucleation isonly possible for temperatures T such that S ( T ) (cid:29) S , and hence a decrease (increase) in the nucle-ation rate. The nucleation rate is strongly dependent on the su-persaturation ratio and hence the temperature; therefore, localover-densities in the atmosphere as a result of an internal grav-ity wave can decrease the local nucleation rate; whereas, localunder-densities can increase the nucleation rate.If there is an established particulate onto which material fromthe atmospheric gas can be absorbed, the growth rate can be alsobe a ff ected by the passage of an internal gravity wave. In this sce-nario, an increase (decrease) in the local density, then the numberof gas particles passing through a target area per unit time alsoincreases (decreases), and so does the number of interactions perunit time that occurs between the seed particulates and the atmo-spheric species.To quantify the impact of the waves on dust nucleation andgrowth, we ran simulations for a range of driving frequenciesand perturbation amplitudes. The density, nucleation and growthresponses are presented in Fig. 5. To obtain results compara-ble across varying frequencies, the values used for plotting weremeasured for τ = T ω , where T ω is the period of the driving oscil-lations, and normalised to their values at static equilibrium in thecentre of the numerical domain. In order to obtain results com-parable across a range of frequencies, the measurements weretaken as the maximum values for ρ , J ∗ , , and γ along a verticalslice of the numerical domain, located at a distance X ( ω ) fromthe centre of the domain so that X ( ω ) = X ref ω ref /ω , where X ref isthe location of the slice picked for a reference case at ω ref .As the perturbation amplitude increases the density ampli-tude of the resulting wave response linearly increases by up to afactor 1 .
5, consistent with the linear regime assumed (top panel,Fig. 5). In contrast, the resulting nucleation rate (middle panel,Fig. 5) increases non-linearly with the perturbation amplitudeby up to a factor 20 in the most favourable scenario ( ω = . N , ρ A = . ρ ). This is a result of the complex non-linear depen-dence of the nucleation rate (Eq. (26)), and not as a consequenceof a non-linear evolution of the internal gravity wave, since oursimulations are conducted in the linear regime only. This impliesthat if the internal gravity wave were to evolve non-linearly, thecorresponding nucleation rate could exhibit an enhanced non-linear response, giving greater nucleation rate values. A similarassertion can be made regarding the mantle growth rate (bot-tom panel, Fig. 5): the mantle growth rate increases by up toa factor 1 . γ ∝ ρ T / = ρ ( γ a − / ; however, the non-linear dependence isless pronounced than that for the nucleation rate. The normalisedgrowth rate γ /γ is independent of the fraction of the surround-ing gas composed of the accreting species f s . Varying the spatialscale of the initial density perturbation does not a ff ect the ampli-tude of the wave response. When ω > N , internal gravity waves cannot propagate sincethe system cannot respond quick enough to the imposed drivenperturbation and the waves are evanescent. When ω ≤ N , in-ternal gravity waves can freely propagate, where the interplaybetween the driving frequency and the local buoyancy frequencyresults in greater response amplitudes for lower driving frequen-cies than for frequencies approaching the natural buoyancy fre-quency of the system (see Fig. 5). As the generated wave prop-agates away from the oscillation source, it encounters regionsof di ff ering background density and hence local buoyancy fre-quency. In response, the wave amplitude, speed, and wavelengthchange in harmony to conserve wave energy. For example, if thewave propagates into regions of lower-density (higher-density),the amplitude of the wave increases (decreases), the wave speedincreases (decreases), and the wavelength decreases (increases)in sympathy. If the wave encounters a region where its frequencyis greater than the local buoyancy frequency, the wave ceases topropagate.
10 5 0 5 10 x /km y / k m / x /km y / k m J *,1 / J *,0 x /km y / k m / Fig. 4.
2D map of ρ /ρ (top panel), J ∗ , / J ∗ , (middle panel) and γ /γ (bottom panel) after one oscillation period. For this example, waves aredriven at ω = . N , and the amplitude of the driving perturbation ρ A = . ρ . The vertical line shows the location of the slice used tomeasure the wave response, at a distance X ( ω ) from the centre of thedomain so that X ( ω ) = X ref ω ref /ω .Article number, page 6 of 10. Parent et al.: Sub-stellar internal gravity waves / = 0.25 × N = 0.50 × N = 0.75 × N = 1.00 × N J * , / J * , A / / Fig. 5.
Plots of ρ /ρ (top panel), J ∗ , / J ∗ , (middle panel), and γ /γ (bottom panel) after one period of an internal gravity wave, as a functionof the amplitude of the density perturbation used to drive the waves.The nucleation plot shows a strong non-linear increase in response toincreased perturbation amplitude, more pronounced with low drivingfrequencies. The impact on dust growth is much weaker, and growthincreases with the driving frequency. The measurements are taken asthe maximum values along a vertical slice of the numerical domain,located at a distance X ( ω ) from the centre of the domain so that X ( ω ) = X ref ω ref /ω . The numerical simulations presented in Fig. 5 are normalisedto equilibrium reference values, aiding in their generalisation toother brown dwarf and exoplanet models. For example, the den-sity response (top plot, Fig 5) and the dust growth rate response(bottom plot, Fig. 5) are indicative of the typical response ex-pected in atmospheric models beyond the exemplar presented.In contrast, the nucleation rate response is more complex andis driven by three main parameters: the ambient gas tempera-ture T , the total gas pressure p , and the number density n x ofthe nucleating species x . These interdependent parameters de-pend upon a number factors, including the chemical composi- tion of the atmosphere and the chemical processes involved (eg.see Helling et al. 2017; Lee et al. 2015, 2018). To investigatethe nucleation rate response, we explored the parameter space ofthe key atmospheric variables, in order to contextualise the re-sults beyond the exemplar atmospheric model considered. Vary-ing T e ff and log g explicitly instead can obfuscate the physicalpicture leading to the underlying cause of the nucleation ratevariations and can be misleading when dealing with microphys-ical processes. Figure 6 shows the contour of the change in nu-cleation rate J ∗ , (normalised to the total nucleation rate J ∗ ) fora range of background temperatures T ∈ [500 K , p ∈ [10 − bar , − bar], and n TiO / n ∈ [10 − , − ] tocover a wide range of values from contemporary models (seeFig. 4, Lee et al. (2015); Fig. 4, Helling et al. (2008a); Fig. 2Stark et al. (2013)).For the top plot of Fig. 6, we set p ≈ × − bar (takenfrom the centre of the domain shown in Fig. 4) and vary n TiO / n and T . For the bottom plot, we held n TiO / n = − (see Fig. 2).We computed J ∗ , using perturbed values of temperature T andpressure p . We obtained those values using the adiabatic equa-tions of state, assuming a density perturbation resulting from thepassage of an internal gravity wave in the most favourable sce-nario ρ /ρ = ± .
48 as shown in Fig. 5. We note that as a re-sult of the wide exploration of parameter space, not all pointsin Fig. 6 correspond to a self-consistently calculated pT atmo-spheric equilibrium state.Figure 6 shows that the increase in nucleation rate J ∗ , makesup a larger portion of the total nucleation rate at higher temper-atures, lower background pressures (top plot), and lower n TiO / n (bottom plot). These conditions are unfavourable to e ffi cient nu-cleation, resulting in very small values of J ∗ , . However, thepassage of an interval gravity wave can produce temporary, lo-calised conditions allowing nucleation to occur at an enhancedrate. While that increase in nucleation may not be large in ab-solute terms, it is much larger than the background values, andleads to a large relative increase. This is consistent with resultsobtained by Helling et al. (2001) for simulated sound waves.To elucidate this point further, we consider a slice of constant n TiO / n in the top plot of Fig. 6. As the temperature increases theequilibrium nucleation rate decreases which inhibits the growthof TiO clusters since nucleation favours cooler temperatures.As a result, the nucleation rate enhancement J ∗ , due to the pas-sage of an internal gravity wave relative to the equilibrium J ∗ , isdiminished leading to an increase in J ∗ , / J ∗ . The opposite occursif the temperature decreases, leading to an increase in J ∗ , / J ∗ . Asimilar response is evident when considering a slice of constant p in the bottom plot of Fig. 6.Further to this, we consider a slice of constant backgroundtemperature in the top plot of Fig. 6. As n TiO / n increases, theequilibrium nucleation rate increases due to the increased su-persaturation ratio S , since the gas-phase TiO molecules aremore likely to cluster and nucleate. As a result, J ∗ , is diminishedrelative to the equilibrium J ∗ , , leading to a decrease in J ∗ , / J ∗ .Similarly, consider a constant slice of temperature T in the bot-tom plot of Fig. 6: increasing the background pressure leads toa greater absolute TiO density n TiO , resulting in an increasedbackground nucleation rate and lower relative increase due to thepassage of an internal gravity wave J ∗ , / J ∗ . If n TiO / n decreasesthe opposite case occurs, leading to an increase in J ∗ , / J ∗ .To place the relative nucleation rate J ∗ , / J ∗ in a wider con-text, in the bottom plot of Fig. 6 we overplot the pT profiles forthe brown dwarf model used for the numerical simulations andof typical and low-gravity L and T dwarfs (see Fig. 1 and Ta-ble 2). Furthermore, in the top plot of Fig. 6 we plot the line Article number, page 7 of 10 & A proofs: manuscript no. parent_2019 of constant temperature corresponding to an atmospheric pres-sure of p ≈ − bar for each of the model atmospheres consid-ered. We note that each model profile has self-consistently cal-culated nucleation rate as a function of atmospheric pressure p ,temperature T , n TiO and atmospheric chemistry that may notcorrespond to a singular point on the parameter space contours.Therefore, these profiles give a helpful indication of the impactof internal gravity waves on the nucleation rate for the varietyof sub-stellar objects considered. For example, we consider thesolid line ( T e ff = , log g = .
0) in the top plot of Fig. 6.For n TiO / n = − , we can deduce J ∗ , / J ∗ ≈ .
95, equiva-lent to J ∗ , ≈ J ∗ , , which is consistent with the data shown inFig. 2. The same result can be obtained using the solid line inthe bottom plot, taking p ≈ × − . The intersection of theconstant pressure line and the model line yields J ∗ , / J ∗ ≈ . J ∗ , ≈ J ∗ , ).In general, Fig. 6 demonstrates that the impact of internalgravity waves on the nucleation rate is significant across the sub-stellar objects considered. The profiles in Fig. 6 show that thestrongest relative increases in nucleation are obtained when theconditions at equilibrium are less suited for e ffi cient nucleation,which leads to any increase caused by the passage of a gravitywave to be comparatively large. This is visible on the top plot,where the warmer L dwarf models (LD: T e ff = , log g = .
5; LGLD: T e ff = , log g = .
0) are linked to strongincrease in nucleation for a wider range of n TiO / n than thecooler T dwarf models (TD: T e ff = , log g = .
5; LGTD: T e ff = , log g = . pT profiles on the bot-tom plots show that for a constant n TiO / n , the warmer L dwar-willf models exhibit a stronger relative increase in nucleationthan T dwarf models.
5. Discussion
This paper has investigated and characterised the e ff ect of linearinternal gravity waves on the evolution of dust clouds in sub-stellar (brown dwarf and gas giant exoplanetary) atmospheresfor the first time. We have shown that in numerical fluid simula-tions, the passage of an internal gravity wave leads to an increaseof dust nucleation by up to a factor 20, and an increase of dustmantle growth rate by up to a factor 1 .
6. Through an explorationof the wider sub-stellar parameter space, we have shown that,in absolute terms, the increase in dust nucleation due to internalgravity waves is stronger in cooler (T dwarfs) and TiO -rich sub-stellar atmospheres. The relative increase, however, is greater inwarm (L dwarf) and TiO -poor atmospheres due to conditionsthat are less suited for e ffi cient nucleation at equilibrium. Thislatter point is important since the stronger the contrast betweenthe perturbed and equilibrium values, the better the chance ofdetecting an observable signal. Recent observations (Maroccoet al. 2014) and models (Hiranaka et al. 2016) suggest that theextreme reddening of some L dwarfs could be due to a dust hazelayer high up in their atmosphere, which could potentially beimpacted by internal gravity waves.The presence of the signature of an internal gravity wave inthe spectra of a brown dwarf could indicate the presence of con-vection deep in the atmosphere or, in the case of a terrestrial ex-oplanetary atmosphere, that the body has a rocky, solid surfacewith relief, since such features are known to trigger the buoyancyoscillation required to generate the waves (Roeten et al. 2019).In such a scenario, the wavelength of the resulting wave couldgive an indication of the scale of the perturbing feature. Fur-ther investigation into the non-linear evolution of internal gravitywaves could potentially yield greater variations in atmospheric n T i O / n . . . . . . . J *,1 / J * contours ( p bar )
500 600 700 800 900 1000 1100 1200 T /K p / b a r . . . . . J *,1 / J * contours ( n TiO / n = 10 ) BDTDLDLGTDLGLD
Fig. 6.
Contour of the change in nucleation rate J ∗ , (normalised to thetotal nucleation rate J ∗ ), as a function of equilibrium temperature T , n TiO / n (top panel) and equilibrium gas pressure p (bottom panel). Theatmospheric profiles for the brown dwarf model used for simulations(BD: T e ff = , log g = . density and nucleation rate. Moreover introducing additional ef-fects, such as the Coriolis e ff ect and dynamical equilibria, andinvestigating their impact on the evolution of internal gravitywaves and the resulting cloud cover, might yield further insightinto inhomogeneous cloud coverage in sub-stellar atmospheres.Additionally, observations of of the photometric variabilityresulting from the propagation of an internal gravity wave couldprovide a novel way of diagnosing the atmospheric gas density.We consider two identical, adjacent vertical atmospheric pro-files. We assume that one profile contains a gas over-density ofamplitude ρ , and a corresponding dust over-density ρ d , that oc-curs over a spatial length scale L , as the result of a propagatinginternal gravity wave. The ratio of the spectral flux density S from both can be expressed in terms of their respective opticaldepths, S S ≈ exp ( τ − τ ) (49) Article number, page 8 of 10. Parent et al.: Sub-stellar internal gravity waves where τ i = (cid:90) D (cid:88) s κ s ρ s d l , (50)and we have assumed that the solid angle subtended by the fea-tures is the same. In a simple approach, we only consider absorp-tion contributions from the gas and the dust in the atmosphere.To simplify matters further, we assume a total mean opacity torepresent the contributions from the gas κ g and from the dust κ d respectively. The non-zero contributions of the optical depth in-tegral over the extent of the atmosphere D , can be approximatedto give, S S ≈ exp [ L ( κ g ρ + κ d ρ d )] . (51)In the linear regime, Figure 5 maps the normalised amplitude ofthe driving density perturbation, ρ A /ρ , to the normalised ampli-tude of the resulting density, ρ /ρ , and nucleation rate, J ∗ / J ∗ ,wave response, ρ ρ = χ ρ A ρ , (52) J ∗ J ∗ = δ ρ A ρ , (53) χ and δ are functions of ρ A /ρ . Therefore, the ratio of flux den-sities can be expressed as, S S ≈ exp [ L ( χκ g ρ A + κ d ρ d )] (54)where ρ d = n d m d , ≈ J ∗ ∆ tm d , ≈ J ∗ δ ρ A ρ ∆ tm d , (55)where m d = π a d ρ m / a d is the radius of a dust grain. To relate the amplitude of the ini-tial density perturbation to the equilibrium atmospheric density,we solve the di ff erential equation for the buoyancy frequency(Eq. (23)) over the length scale of the perturbation L , giving ρ A ≈ ρ exp ( − qN L /g ) , (56)where we have assumed that the N is approximately constantacross L , and qN is the wave frequency with q ∈ (0 ,
1] (weassume q = ρ ≈ exp ( N L /g ) L κ g χ [ln ( S / S ) − L κ d ρ d ] , (57)where ρ d ≈ J ∗ δ ∆ tm d exp ( − N L /g ) . (58)This expression allows us to estimate the density of a sub-stellar atmosphere based on potential observations of S / S ona timescale consistent with the buoyancy frequency. As an ex-ample, for the order of magnitude values list in Table 3, Eq. (58)gives an estimation of ρ ≈ − kg m − , which is consistent with contemporary atmospheric numerical models. This demon-strates that from observations of S / S , L , and the timescaleof variation, an estimation of the atmospheric density of a sub-stellar atmosphere can be made. A more in-depth analysis couldinvolve calculating synthetic spectra showing the impact of thegravity wave that could be expected from observations, and willbe considered in a further paper. Acknowledgements.
The authors are grateful to the anonymous referee for con-structive comments and suggestions that have improved this paper. A.P. is grate-ful for funding and support received from Abertay University as part of theRLINCS studentship programme. C.R.S. is grateful for funding from the RoyalSociety via grant number RG160840 and from the Carnegie Trust for the Univer-sities of Scotland via research incentive grant number RIG007788. G.K.H. Leeacknowledges support from the University of Oxford and CSH Bern throughthe Bernoulli fellowship and support from the European community through theERC advanced grant project EXOCONDENSE (PI: R.T. Pierrehumbert).
References
Alvan, L., Brun, A. S., & Mathis, S. 2014, A&A, 565, A42Biller, B. 2017, Astronomical Review, 13, 1Buenzli, E., Apai, D., Radigan, J., Reid, I. N., & Flateau, D. 2014, ApJ, 782, 77Freytag, B., Allard, F., Ludwig, H. G., Homeier, D., & Ste ff en, M. 2010, A&A,513, A19Gail, H.-P., Keller, R., & Sedlmayr, E. 1984, A&A, 133, 320Garaud, P., Medrano, M., Brown, J. M., Mankovich, C., & Moore, K. 2015, ApJ,808Gillon, M., Triaud, A. H. M. J., Jehin, E., et al. 2013, A&A, 555, L5Hauschildt, P. H. & Baron, E. 1999, Journal of Computational and Applied Math-ematics, 109, 41Helling, C., Ackerman, A., Allard, F., et al. 2008a, MNRAS, 391, 1854Helling, C. & Casewell, S. 2014, A&Ar, 22, 80Helling, C., Klein, R., Woitke, P., Nowak, U., & Sedlmayr, E. 2004, A&A, 423,657Helling, C., Oevermann, M., Lüttke, M. J. H., Klein, R., & Sedlmayr, E. 2001,A&A, 376, 194Helling, C., Tootill, D., Woitke, P., & Lee, G. 2017, A&A, 603, A123Helling, C. & Woitke, P. 2006, A&A, 455, 325Helling, C., Woitke, P., & Thi, W. F. 2008b, A&A, 485, 547Hiranaka, K., Cruz, K. L., Douglas, S. T., Marley, M. S., & Baldassare, V. F.2016, ApJ, 830, 96Lee, G., Dobbs-Dixon, I., Helling, C., Bognar, K., & Woitke, P. 2016, A&A, 594,A48Lee, G., Helling, C., Giles, H., & Bromley, S. T. 2015, A&A, 575, A11Lee, G. K. H., Blecic, J., & Helling, C. 2018, A&A, 614, A126Marley, M. S. & Robinson, T. D. 2015, ARA&A, 53, 279Marley, M. S., Saumon, D., & Goldblatt, C. 2010, ApJ, 723, L117Marocco, F., Day-Jones, A. C., Lucas, P. W., et al. 2014, MNRAS, 439, 372Mittal, S. 2014, International journal of high performance computing and net-working, 7, 292Robinson, T. D. & Marley, M. S. 2014, ApJ, 785, 158Rodríguez-Barrera, M. I., Helling, C., & Wood, K. 2018, A&A, 618, A107Roeten, K. J., Bougher, S. W., Benna, M., et al. 2019, Journal of GeophysicalResearch (Planets), 124, 3283Saumon, D. & Marley, M. S. 2008, ApJ, 689, 1327Stark, C. R. & Diver, D. A. 2018, A&A, 611, A91Stark, C. R., Helling, C., & Diver, D. A. 2015, A&A, 579, A41Stark, C. R., Helling, C., Diver, D. A., & Rimmer, P. B. 2013, ApJ, 776, 11Sutherland, B. R. 2010, Internal Gravity Waves (Cambridge Univerisity Press)Tremblin, P., Amundsen, D. S., Chabrier, G., et al. 2016, ApJ, 817, L19Tremblin, P., Amundsen, D. S., Mourier, P., et al. 2015, ApJ, 804, L17Tremblin, P., Padioleau, T., Phillips, M. W., et al. 2019, ApJ, 876, 144Vallis, G. K. 2017, Atmospheric and Oceanic Fluid Dynamics (Cambridge Uni-versity Press)Vetterling, W. T., Teukolsky, S. A., Press, W. H., & Flannery, B. P. 1992, Numer-ical recipes in C, 2nd edn. (Cambridge University Press)Vos, J. M., Allers, K., Apai, D., et al. 2019, arXiv e-prints, arXiv:1903.06691Witte, S., Helling, C., Barman, T., Heidrich, N., & Hauschildt, P. 2011, A&A,529, A44Witte, S., Helling, C., & Hauschildt, P. 2009, A&A, 506, 1367 Article number, page 9 of 10 & A proofs: manuscript no. parent_2019
Table 3.
Example order of magnitude values used in Eq. (58)
Parameter Value Notes L m See Freytag et al. (2010); Marley & Robinson (2015) κ g − m kg See Fig. 13, Lee et al. (2016) for λ ≈ µ m κ d − m kg See Fig. 13, Lee et al. (2016) for λ ≈ µ m ∆ t = T N =
10 s See Fig. 3 χ ρ A /ρ = . ω/ N = . δ
30 Fig. 5 for ρ A /ρ = . ω/ N = . J ∗ m − s − Eq. (26) at equilibrium, see also Fig. 1 in Helling et al. (2008b) S / S .
001 see Figs 7 −
10 in Buenzli et al. (2014) m d − kg for a d = − m, ρ = kg m − g m s − log g = N = π / ∆ tt