Note

Download full notebook here

Chemistry: Young Planet Spectroscopy

What you will learn:

  1. What happens to the spectra of planet atmospheres as they cool?

  2. What is happening from a chemical standpoint to cause these spectral features?

  3. What molecules are good temperature probes?

What you should know:

  1. What do formation models predict for the effective temperatures of young planets across different masses?

  2. Given identical luminosity and age, can formation scenarios and mass be determined?

  3. How do we dissect spectroscopy of planet atmospheres in order to infer atmospheric physical properties such as abundance and climate profiles?

[1]:

import warnings
warnings.filterwarnings(action='ignore')
import picaso.justdoit as jdi
import picaso.justplotit as jpi
jpi.output_notebook()
import os
import pandas as pd
import numpy as np
#point to your sonora profile grid that you untared (see above cell #2)
sonora_profile_db = '/data/sonora_profile/'
Loading BokehJS ... 
[2]:

wave_range = [0.6,5] #don't worry we will play around with this more later
opa = jdi.opannection(wave_range=wave_range)

What happens to an observable atmospheric spectrum as the planet ages

In the previous workbook, we learned how to analyze spectra. Now we will compare a sequence of spectra as a function of age in order to gain an intuition for major transitions expected as a function of age (and by proxy, temperature and gravity).

[3]:

case_study = jdi.evolution_track(mass=8, age='all')

There are a lot of grid points on this evolutionary track (150!). Let’s pick the “hot” start case as it offers a more dramatic cooling event. This will enable us to learn about the chemical transitions that happen as a planet cools.

[4]:

case_study['hot']

[4]:
age_years Teff grav_cgs logL R_cm
1 1.006340e+06 2276.0 3271.9 30.773 1.760390e+10
2 1.012680e+06 2274.6 3287.2 30.770 1.756300e+10
3 1.016060e+06 2273.9 3295.3 30.768 1.754130e+10
4 1.020340e+06 2272.9 3305.5 30.766 1.751410e+10
5 1.025760e+06 2271.7 3318.5 30.763 1.747990e+10
... ... ... ... ... ...
146 1.567600e+10 200.5 19931.3 25.768 7.132490e+09
147 1.657590e+10 197.2 19985.4 25.738 7.122830e+09
148 1.777590e+10 193.2 20052.3 25.701 7.110940e+09
149 1.937580e+10 188.3 20133.4 25.655 7.096600e+09
150 2.000100e+10 186.6 20162.9 25.638 7.091400e+09

150 rows × 5 columns

Let’s take a feasible subset of these. I will choose ten, though if you are curious, or want to do the full track go for it!

[5]:

i_to_compute = case_study['hot'].index[0::15]#take every 15th value

Let’s run PICASO in a loop with the different effective temperatures and gravities

[6]:

yph = jdi.inputs()
#let's keep the star fixed
yph.star(opa, 5000,0,4.0,radius=1, radius_unit=jdi.u.Unit('R_sun'))
yph.phase_angle(0)

#Let's stick the loop in right here!
hot_output={} #easy mechanism to save all the output
for i in i_to_compute:
    Teff = case_study['hot'].loc[i,'Teff']
    grav = case_study['hot'].loc[i,'grav_cgs']
    yph.gravity(gravity= grav,
                gravity_unit=jdi.u.Unit('cm/s**2'))
    yph.sonora(sonora_profile_db,  Teff)
    hot_case = yph.spectrum(opa,calculation='thermal', full_output=True)
    hot_output[f'{Teff}_{grav}'] = hot_case

Let’s plot the sequence!!

[7]:

wno,spec=[],[]
fig = jpi.figure(height=500,width=600, y_axis_type='log',
                 x_axis_label='Wavelength(um)',y_axis_label='Flux (erg/s/cm2/cm)')
for i,ikey in enumerate(hot_output.keys()):
    x,y = jdi.mean_regrid(hot_output[ikey]['wavenumber'],
                          hot_output[ikey]['thermal'], R=150)

    t,g=tuple(ikey.split('_'));g=int(np.log10(float(g))*1000)/1000
    a=fig.line(1e4/x,y,color=jpi.pals.Spectral11[::-1][i],line_width=3,
               legend_label=f'Teff={t}K,logg={g}')
fig.legend.location='bottom_right'

jpi.show(fig)

There is rich information encoded in these spectra. In order to fully grasp what is going on, it is important to understand the chemistry.

What molecules are most important to planetary spectroscopy?

In the previous exercise we focused on look at the specific molecular contributions of two distinct cases. Therefore, we were focused on abundances as a function of pressure. Here we want you to get a handle on bulk abundance properties as a function of effective temperature. So we are going to collapse the pressure axis by taking the “median” value of each abundance array. By doing so, we want to see what the ~10 most abundant molecules are in each of these 10 spectra.

[8]:

#remember the mixing ratios (or abundances) exist in this pandas dataframe
hot_output[ikey]['full_output']['layer']['mixingratios'].head()
#but this is too many molecules to keep track of for every single spectrum

[8]:
e- H2 H H+ H- H2- H2+ H3+ He H2O ... SiO MgH OCS Li LiOH LiH LiCl graphite Al CaH
0 4.505501e-38 0.836890 4.505501e-38 4.505501e-38 4.505501e-38 4.505501e-38 4.505501e-38 4.505501e-38 0.16265 1.992481e-25 ... 4.505501e-38 4.505501e-38 4.505501e-38 4.505501e-38 4.505501e-38 4.505501e-38 4.505501e-38 4.505501e-38 5.325528e-57 5.325528e-57
1 4.505500e-38 0.836890 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 0.16265 4.362606e-25 ... 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.553439e-57 4.553439e-57
2 4.505500e-38 0.836889 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 0.16265 1.160365e-24 ... 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 3.907596e-57 3.907596e-57
3 4.505500e-38 0.836889 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 0.16265 3.579045e-24 ... 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 3.358949e-57 3.358949e-57
4 4.505500e-38 0.836888 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 0.16265 1.400859e-23 ... 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 4.505500e-38 2.895678e-57 2.895678e-57

5 rows × 40 columns

[9]:

relevant_molecules=[]
for i,ikey in enumerate(hot_output.keys()):
    abundances = hot_output[ikey]['full_output']['layer']['mixingratios']

    #first let's get the top 10 most abundance species in each model bundle we ran
    median_top_10 = abundances.median().sort_values(ascending=False)[0:10]
    relevant_molecules += list(median_top_10.keys())

#taking the unique of  relevant_molecules will give us the molecules we want to track
relevant_molecules = np.unique(relevant_molecules)

print(relevant_molecules)

['Al' 'C2H6' 'CH4' 'CO' 'CO2' 'Fe' 'H' 'H2' 'H2O' 'H2S' 'He' 'K' 'LiCl'
 'N2' 'NH3' 'Na' 'PH3' 'SiO']

Now that we have condensed this to a meaningful set of molecules, we can proceed to plot the sequence

Side note: You might try to see if the technique of taking the “median” yields the same results as “max” or “mean”. This gives some insight into how dynamic moleculare abundances are as a function of pressure

Where in temperature space do chemical transitions seem to take place?

[10]:

fig = jpi.figure(height=500,width=700, y_axis_type='log',
                 y_range=[1e-15,1],x_range=[200,2600],
                 x_axis_label='Planet Effective Temperature',y_axis_label='Abundance')

#now let's go back through our models and plot the abundances as a function of teff
relevant_molecules={i:[] for i in relevant_molecules}
for i,ikey in enumerate(hot_output.keys()):
    abundances = hot_output[ikey]['full_output']['layer']['mixingratios'].median()

    #save each abundance
    for i in relevant_molecules.keys():
        relevant_molecules[i] += [abundances[i]]

#last loop to plot each line
for i,ikey in enumerate( relevant_molecules.keys()):
    fig.line(case_study['hot'].loc[i_to_compute,'Teff'], relevant_molecules[ikey],
               color=jpi.pals.Category20[20][i],line_width=3,legend_label=ikey)
fig.legend.location='bottom_right'
jpi.show(fig)

There is a lot happening but let’s break it down in very broad digestible categories. I will ask you to look back to the spectra that you made in the first tutorial. However in some cases, those spectra might not be computed at high enough effective temperatures to explore the molecular contribution. In those cases, use the techniques you learned from the previous exercise (jdi.get_contribution) to answer the questions below:

Universally abundant molecules:

  • Which are the few highest abundance molecules/elements that exist across all temperature?

  • In what ways do these molecules/elements contribute to planetary spectra?

Carbon-bearing species (CO2, CH4, CO, C2H6):

  • Which molecules are good temperature indicators, meaning they only exist in certain temperature regimes?

  • For the molecules that are good temperature indicators, where do their transitions occur? Keep these numbers archived in the back of your brain as they are great to have for intuition

  • Do these molecules make a significant contribution to the spectra? Are they hard or easy to detect? At what wavelengths?

Besides Carbon, what other non-metal-based molecules are dominant?

  • Are any of them indicators of high or low temperature?

  • Do any of them exhibit an interplay that is similar to that of the CH4/CO transition?

  • Do these molecules make a significant contribution to the spectra? Are they hard or easy to detect? At what wavelengths?

What Alkali-based molecules/elements are dominant?

  • At what temperatures do these molecules/elements begin to appear?

  • Do these molecules make a significant contribution to the spectra? Are they hard or easy to detect? At what wavelengths?

What Metal-based species are dominant?

  • At what temperatures do these molecules/elements begin to appear?

  • Do these molecules make a significant contribution to the spectra? Are they hard or easy to detect? At what wavelengths?

SYNTHESIZE:

  • Across all these molecules, what are the few most critical temperature transitions? According to our two formation scenarios, what planet age does these effective temperatures correspond to?