One-Dimensional Climate Models: Brown Dwarfs w/ Clouds in Chemical Disequilibrium
In this tutorial you will learn how to run 1-D climate models for brown dwarfs with Virga clouds and chemical disequilibrium. For a more in depth look at the climate-cloud code check out Mang et al. 2026 (note this should also be cited if using this code/tutorial).
You should already be familiar with running 1-D climate models with running a simple clear brown dwarf model and running clouds in equilibrium
What you need to have downloaded for clouds to work:
[1]:
from virga import justplotit as vpi
from virga import justdoit as vdi
from picaso import justplotit as jpi
from picaso import justdoit as jdi
import os
import warnings
warnings.filterwarnings('ignore')
jpi.output_notebook()
import astropy.units as u
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
from bokeh.plotting import show, figure
import xarray as xr
import pickle
WARNING: Failed to load Vega spectrum from /data/reference_data/picaso/ref4/stellar_grids/calspec/alpha_lyr_stis_011.fits; Functionality involving Vega will be severely limited: FileNotFoundError(2, 'No such file or directory') [stsynphot.spectrum]
[2]:
# #sonora bobcat cloud free structures file
# sonora_profile_db = '/data/sonora_bobcat/structure/structures_m+0.0'
# sonora_profile_db = '/Users/nbatalh1/Documents/data/sonora_bobcat/structures_m+0.0'
sonora_profile_db = os.path.join(os.getenv('picaso_refdata'),'sonora_grids','bobcat')
We are going to initialize our climate run like other disequilibrium runs without using the pre-weighted ck tables
[3]:
cl_run = jdi.inputs(calculation="browndwarf", climate = True) # start a calculation
#note you need to put the climate keyword to be True in order to do so
# now you need to add these parameters to your calculation
teff= 400 # Effective Temperature of your Brown Dwarf in K
grav = 1000 # Gravity of your brown dwarf in m/s/s
cl_run.gravity(gravity=grav, gravity_unit=u.Unit('m/(s**2)')) # input gravity
cl_run.effective_temp(teff) # input effective temperature
gases_fly = ['CO','CH4','H2O','NH3','CO2','N2','HCN','H2','PH3','C2H2','Na','K','TiO','VO','FeH']
opacity_ck = jdi.opannection(method='resortrebin',preload_gases=gases_fly)# grab your opacities
[4]:
nlevel = 91 # number of plane-parallel levels in your code
# Here we're going to start with a cloudfree Sonora Elf Owl model
pressure,temp_guess = np.loadtxt("profilegrid_kz_1d9_qt_onfly_400_grav_1000_mh_+0.0_cto_1.0.dat",
usecols=[1,2],unpack=True, skiprows = 1)
[5]:
rcb_guess = 45 # top most level of guessed convective zone
# for the sake of time of this tutorial, I set it to 40 because I know where it should be in this case. In general for clouds it is better
# to start deeper in the atmosphere and work your way up. It just takes more time.
# Here are some other parameters needed for the code.
rfacv = 0.0 #we are focused on a brown dwarf so let's keep this as is
# this is the Kzz parameter
kzval = pressure*0+1e9
[6]:
virga_dir = os.path.join(os.getenv('picaso_refdata'),'virga') # path to virga directory with the Mie scattering files
All requirements needed for this notebook will be met from properly following the previous notebook for clouds in equilibrium chemistry
New PICASO code parameters: We’ve introduced other optional chemical treatments into the disequilibrium runs that should be included in the atmosphere() function
vol_rainout: (True/False) Follows rainout chemistry (ie. Sonora Bobcat) even in disequilibrium. Default = Falsecold_trap: (True/False) If True, the abundances of volatile species like H2O, NH3, and CH4 will not be allowed to increase after they begin to rainout. Default = Falseno_ph3: (True/False) If True, completely remove PH3 from the atmosphere. Default = False
These parameters are included to increase flexibility for you when you want to generate your own model! To play around with these different parameters you can look at the chemistry notebooks
Note that all disequilibrium models default to follow the kinetic CO2 prescription described in Zahnle & Marley 2014 where the abundance quenches with respect to CO in disequilibrium. This is a correction applied to Sonora Elf Owl v2
Our recommendation to follow for the latest release of the Sonora family of models are as follows:
vol_rainout= Truecold_trap= Trueno_ph3= True
[7]:
cl_run.inputs_climate(temp_guess= temp_guess, pressure= pressure,
rcb_guess=rcb_guess, rfacv = rfacv)
cl_run.atmosphere(mh=1, cto_relative=1, chem_method='visscher',
quench=True, cold_trap = True, no_ph3 = True, vol_rainout=True)#
cl_run.inputs['atmosphere']['profile']['kz']=kzval
# Let's add in a few different cloud species that we expect to condense in the deeper parts of the atmosphere of this 400 K brown dwarf
df_cld_guess=cl_run.virga(condensates=['Cr', 'MnS', 'Na2S', 'ZnS', 'KCl'], fsed=2.0, directory=virga_dir, latent_heat=True)
In this case, with a good starting radiative-convective boundary guess and a good starting PT profile, this run should take ~10 minutes to compute. If it takes significantly longer longer there might be something wrong. Most times though, cloudy runs do take much longer on timescales of 30 minutes to an hour depending on how deep your intial radiaitve-convective zone boundary is.
[8]:
out = cl_run.climate(opacity_ck, save_all_profiles = True, with_spec=True,
diseq_chem = True, self_consistent_kzz =False)
SUMMARY
-------
Clouds: True
quench True
cold_trap True
vol_rainout True
no_ph3 True
Moist Adiabat: False
Kzz for chem: [1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09
1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09
1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09
1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09
1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09
1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09
1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09
1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09
1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09 1.e+09
1.e+09]
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
I dont have fluxes, let me compute them
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(67), 'NH3-N2': np.int64(75), 'HCN': np.int64(74), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 2.3383512429196576 Tau tolerance is 0.4676702485839315
Iteration number 0 , min , max temp 96.96845780821303 3475.321434279068 , flux balance -0.04740854436036882
Iteration number 1 , min , max temp 101.24659023080862 3475.36894254416 , flux balance -0.04735523150179005
Iteration number 2 , min , max temp 105.03087116078686 3475.4275795100552 , flux balance -0.047289422623226114
Iteration number 3 , min , max temp 108.60553636880357 3475.499806729022 , flux balance -0.04720834975209569
Iteration number 4 , min , max temp 111.99961628019263 3475.588417198293 , flux balance -0.04710886987119428
Iteration number 5 , min , max temp 115.24206355042318 3475.696588257649 , flux balance -0.04698740421452066
Iteration number 6 , min , max temp 118.35729402850885 3475.8279174385566 , flux balance -0.046839896210493524
Iteration number 7 , min , max temp 121.36464778807523 3475.9864416863534 , flux balance -0.0466617874889959
Iteration number 8 , min , max temp 124.23523050878141 3476.176639689295 , flux balance -0.04644801213706861
Iteration number 9 , min , max temp 126.99464451169433 3476.4034161062636 , flux balance -0.04619301045341567
Iterations exceeded max_inner_iterations ! sorry
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(67), 'NH3-N2': np.int64(75), 'HCN': np.int64(74), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 2.894833262333518 Tau tolerance is 1.5002362458556866
Big iteration is 126.99464451169433 0
Iteration number 0 , min , max temp 129.32718177825694 3477.354369876842 , flux balance -0.06831280324827038
Iteration number 1 , min , max temp 131.61390420466842 3478.375648763716 , flux balance -0.06704290147344737
Iteration number 2 , min , max temp 133.85706496975297 3479.471627121645 , flux balance -0.06567645701785789
Iteration number 3 , min , max temp 136.05903308070995 3480.647519443145 , flux balance -0.06420641534724399
Iteration number 4 , min , max temp 138.22217326604502 3481.9091217756595 , flux balance -0.06262561558066729
Iteration number 5 , min , max temp 140.34874467939682 3483.261573249669 , flux balance -0.06092682317481027
Iteration number 6 , min , max temp 142.44085756388617 3484.7099450262185 , flux balance -0.05910276349288195
Iteration number 7 , min , max temp 144.50044413070765 3486.2592067163882 , flux balance -0.057146152864755176
Iteration number 8 , min , max temp 146.5292639746011 3487.914193694075 , flux balance -0.05504972862672992
Iteration number 9 , min , max temp 148.52890281384575 3489.6795763234772 , flux balance -0.05280627697286986
Iterations exceeded max_inner_iterations ! sorry
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(67), 'NH3-N2': np.int64(75), 'HCN': np.int64(74), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 1.6901220477022516 Tau tolerance is 2.347750271661543
Big iteration is 148.52890281384575 1
Iteration number 0 , min , max temp 149.23849221711592 3492.241967859119 , flux balance -0.06807666298454727
Iteration number 1 , min , max temp 149.98283496386028 3494.9483906309274 , flux balance -0.06445932601716989
Iteration number 2 , min , max temp 150.76196512264426 3497.8002961188936 , flux balance -0.06062538352969672
Iteration number 3 , min , max temp 151.57568742501178 3500.8017607198994 , flux balance -0.05656698102029992
Iteration number 4 , min , max temp 152.4235936385299 3503.956232285873 , flux balance -0.05227649098758673
Iteration number 5 , min , max temp 153.30508485472026 3507.2640252143005 , flux balance -0.04774664725967454
Iteration number 6 , min , max temp 154.2193921693227 3510.7282077453347 , flux balance -0.04297059190021088
Iteration number 7 , min , max temp 155.16560220651874 3514.35231263811 , flux balance -0.03794186622282363
Iteration number 8 , min , max temp 156.14269306750543 3518.1372688919046 , flux balance -0.03265451049620892
Iteration number 9 , min , max temp 157.14958095376866 3522.0781230399216 , flux balance -0.027103229461770094
Iterations exceeded max_inner_iterations ! sorry
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(74), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 1.629480519820171 Tau tolerance is 2.8622636286699366
Big iteration is 157.14958095376866 2
Iteration number 0 , min , max temp 158.20520176062965 3528.7223552814166 , flux balance -0.044322519490276324
Iteration number 1 , min , max temp 159.28416824686136 3535.5412936180155 , flux balance -0.03365639918913065
Iteration number 2 , min , max temp 160.3857729477442 3542.5429116836694 , flux balance -0.022567072576480662
Iteration number 3 , min , max temp 161.51007372157017 3549.694408644163 , flux balance -0.011060980961848295
Iteration number 4 , min , max temp 162.56894874329117 3556.3892495993673 , flux balance -0.00012847368309654324
Iteration number 5 , min , max temp 162.596767619879 3556.4593833333606 , flux balance -2.149719236143551e-06
In t_start: Converged Solution in iterations 5
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 1.6105277852101203 Tau tolerance is 2.774178065992
Big iteration is 162.596767619879 3
Iteration number 0 , min , max temp 161.72614447237302 3563.473760635926 , flux balance -0.003064162500728795
Iteration number 1 , min , max temp 161.48770025134726 3565.1200935869538 , flux balance -4.662037161922596e-05
In t_start: Converged Solution in iterations 1
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(72), 'CO2': np.int64(66), 'NH3-N2': np.int64(74), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 2.7199872851562725 Tau tolerance is 2.5595869524559705
Big iteration is 161.48770025134726 4
Iteration number 0 , min , max temp 160.77829550401222 3559.19097530382 , flux balance 0.0002587168407406222
Iteration number 1 , min , max temp 160.75234352165708 3559.0684301847386 , flux balance 4.324940175551415e-06
In t_start: Converged Solution in iterations 1
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 1.2365564426934643 Tau tolerance is 2.497291805497537
Profile converged before max_outer_iterations
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
I dont have fluxes, let me compute them
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(74), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 0.4384594831051678 Tau tolerance is 2.4866248646847966
Iteration number 0 , min , max temp 161.62853075098624 3540.749425332699 , flux balance -0.02339009989668019
Iteration number 1 , min , max temp 162.20771912046982 3554.9865617672804 , flux balance -0.00015999543382353573
Iteration number 2 , min , max temp 162.22670514827934 3555.0626825877066 , flux balance -2.696922981646181e-06
In t_start: Converged Solution in iterations 2
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 0.3341489492856802 Tau tolerance is 2.45615391790127
Big iteration is 162.22670514827934 0
Iteration number 0 , min , max temp 162.11878483778162 3562.60870553068 , flux balance -0.0001633659850910556
Iteration number 1 , min , max temp 162.11861013300992 3562.6891257287743 , flux balance -2.761035487855922e-06
In t_start: Converged Solution in iterations 1
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(74), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 2.6958482435308806 Tau tolerance is 2.4592717586001878
Big iteration is 162.11861013300992 1
Iteration number 0 , min , max temp 162.39710180101758 3564.7863209402044 , flux balance -6.415946407050036e-05
Iteration number 1 , min , max temp 162.40810074766654 3564.817908681757 , flux balance -1.0850696267400674e-06
In t_start: Converged Solution in iterations 1
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(72), 'CO2': np.int64(66), 'NH3-N2': np.int64(74), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 2.806590885887198 Tau tolerance is 2.468369215411511
Big iteration is 162.40810074766654 2
Iteration number 0 , min , max temp 161.32581677992974 3555.377397067944 , flux balance 0.0004738827615532695
Iteration number 1 , min , max temp 161.2704117458016 3555.146799486102 , flux balance 8.169418042930362e-06
In t_start: Converged Solution in iterations 1
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 0.604539286133075 Tau tolerance is 2.5112499696099837
Profile converged before max_outer_iterations
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
I dont have fluxes, let me compute them
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(74), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 0.6821493189082082 Tau tolerance is 2.5048706714609708
Iteration number 0 , min , max temp 162.20662407480413 3546.134241251439 , flux balance -0.009996650557662475
Iteration number 1 , min , max temp 162.46098914363918 3551.9681490367116 , flux balance -0.00014264397212528971
In t_start: Converged Solution in iterations 1
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 0.15887140053292104 Tau tolerance is 2.4746914751648603
Big iteration is 162.46098914363918 0
Iteration number 0 , min , max temp 162.3741348811857 3558.1753141198988 , flux balance -0.00016146945626562882
Iteration number 1 , min , max temp 162.3730376039339 3558.255399632255 , flux balance -2.9435604946686962e-06
In t_start: Converged Solution in iterations 1
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 2.78029794839878 Tau tolerance is 2.4921816640268846
Big iteration is 162.3730376039339 1
Iteration number 0 , min , max temp 162.65359703601467 3559.58988898621 , flux balance -4.741708073658078e-05
In t_start: Converged Solution in iterations 0
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 2.75838902967864 Tau tolerance is 2.4942152703857747
Big iteration is 162.65359703601467 2
Iteration number 0 , min , max temp 162.672888674599 3559.9267826254095 , flux balance -1.1782269056142802e-05
In t_start: Converged Solution in iterations 0
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 0.6278858199177151 Tau tolerance is 2.4775780341523688
Profile converged before max_outer_iterations
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
I dont have fluxes, let me compute them
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(74), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 0.6102049569115042 Tau tolerance is 2.489977253409261
Iteration number 0 , min , max temp 162.6222268571807 3553.8674130641857 , flux balance -0.00020310568651806302
Iteration number 1 , min , max temp 162.62400806146886 3553.9720495816305 , flux balance -3.955145947996608e-06
In t_start: Converged Solution in iterations 1
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 0.05819631100569156 Tau tolerance is 2.471156150234398
Big iteration is 162.62400806146886 0
Iteration number 0 , min , max temp 162.77278687707906 3557.1044231064884 , flux balance -0.00010382410496234265
In t_start: Converged Solution in iterations 0
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 2.7303713246789267 Tau tolerance is 2.4715060648744394
Big iteration is 162.77278687707906 1
Iteration number 0 , min , max temp 162.80443420270896 3557.7916300857896 , flux balance -2.6423614894042622e-05
In t_start: Converged Solution in iterations 0
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 2.7049969355695445 Tau tolerance is 2.479214030033425
Big iteration is 162.80443420270896 2
Iteration number 0 , min , max temp 162.8217720856083 3557.8934571515056 , flux balance -4.046074969584501e-06
In t_start: Converged Solution in iterations 0
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 0.58970841753057 Tau tolerance is 2.4766261114065564
Profile converged before max_outer_iterations
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
I dont have fluxes, let me compute them
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 0.5690275963592235 Tau tolerance is 2.4823097783564534
Iteration number 0 , min , max temp 162.66989584887617 3555.770925433346 , flux balance -0.0001746733189915663
In t_start: Converged Solution in iterations 0
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 0.014839203156924441 Tau tolerance is 2.4676985402020115
Big iteration is 162.66989584887617 0
Iteration number 0 , min , max temp 162.89795288780587 3557.11081593958 , flux balance -5.2949047854077814e-05
In t_start: Converged Solution in iterations 0
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 2.7218614261042493 Tau tolerance is 2.4554815817245617
Big iteration is 162.89795288780587 1
Iteration number 0 , min , max temp 162.86544785773782 3557.4211640694766 , flux balance -1.3241178159069385e-05
In t_start: Converged Solution in iterations 0
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 2.7208301512919384 Tau tolerance is 2.458471198021276
Big iteration is 162.86544785773782 2
Iteration number 0 , min , max temp 162.89628609960633 3557.4806980635317 , flux balance -2.6159253041514187e-06
In t_start: Converged Solution in iterations 0
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 0.5688019430860134 Tau tolerance is 2.4830224536221053
Profile converged before max_outer_iterations
final [0, 41, 89, 0, 0, 0]
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
I dont have fluxes, let me compute them
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 0.5537814697981531 Tau tolerance is 2.486043384223722
Iteration number 0 , min , max temp 162.90468791437218 3557.6407360976737 , flux balance -6.644996133798946e-06
In t_start: Converged Solution in iterations 0
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 0.0056910613814649125 Tau tolerance is 2.4728330844909716
Big iteration is 162.90468791437218 0
Iteration number 0 , min , max temp 162.9005488164091 3557.68835437999 , flux balance -2.026080106926797e-06
In t_start: Converged Solution in iterations 0
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 2.721798136478684 Tau tolerance is 2.4401367485287078
Big iteration is 162.9005488164091 1
We are already at a root, tolf , test = 5e-05 , 4.130269885944352e-05
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 2.731300132243704 Tau tolerance is 2.440990221458394
Big iteration is 162.9005488164091 2
We are already at a root, tolf , test = 5e-05 , 4.470313807415311e-05
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
I am updating kzz. This could be either because clouds were requested (which always use self consistent kzz) or because self consistent kzz profiles were requested for chemistry.
Computed quenched levels at {'CO-CH4-H2O': np.int64(73), 'CO2': np.int64(66), 'NH3-N2': np.int64(75), 'HCN': np.int64(73), 'PH3': np.int64(72)}
Quench=True; Adjusting quench chemistry
vol_rainout=True; Adjusting volatile rainout
cold_trap=True; Adjusting cold trap
no_ph3=True; Goodbye PH3!
Doing clouds: Max TAUCLD diff is 0.5587894277899172 Tau tolerance is 2.4921325577304962
Profile converged before max_outer_iterations
YAY ! ENDING WITH CONVERGENCE
Plot the P-T Profile
Now we can plot the results, first let’s grab the condensation curve for all our cloud species
[9]:
kcl_cond_p, kcl_cond_t = vdi.condensation_t('KCl', 1, 2.2, pressure = out['pressure'])
na2s_cond_p, na2s_cond_t = vdi.condensation_t('Na2S', 1, 2.2, pressure = out['pressure'])
mns_cond_p, mns_cond_t = vdi.condensation_t('MnS', 1, 2.2, pressure = out['pressure'])
zns_cond_p, zns_cond_t = vdi.condensation_t('ZnS', 1, 2.2, pressure = out['pressure'])
cr_cond_p, cr_cond_t = vdi.condensation_t('Cr', 1, 2.2, pressure = out['pressure'])
[10]:
pressure_bobcat,temp_bobcat = np.loadtxt(jdi.os.path.join(
sonora_profile_db,f"t{teff}g{grav}nc_m0.0.cmp.gz"),
usecols=[1,2],unpack=True, skiprows = 1)
pressure_elfowl,temp_elfowl = np.loadtxt(f"profilegrid_kz_1d9_qt_onfly_400_grav_1000_mh_+0.0_cto_1.0.dat",
usecols=[1,2],unpack=True, skiprows = 1)
plt.figure(figsize=(10,8))
plt.ylabel("Pressure [Bars]")
plt.xlabel('Temperature [K]')
plt.xlim(0,max(out['temperature'])+50)
plt.ylim(3e3,1e-3)
plt.semilogy(temp_bobcat,pressure_bobcat,color="k",linestyle="--",linewidth=3,label="Sonora Bobcat")
plt.semilogy(temp_elfowl,pressure_elfowl,color='r',linewidth=3, label="Sonora Elf Owl")
plt.semilogy(out['temperature'],out['pressure'],label="Our Cloudy Run")
plt.plot(kcl_cond_t,kcl_cond_p, color ='gray', label = 'KCl Condensation Curve')
plt.plot(na2s_cond_t,na2s_cond_p, color ='orange', label = 'Na2S Condensation Curve')
plt.plot(mns_cond_t,mns_cond_p, color ='purple', label = 'MnS Condensation Curve')
plt.plot(zns_cond_t,zns_cond_p, color ='green', label = 'ZnS Condensation Curve')
plt.plot(cr_cond_t,cr_cond_p, color ='blue', label = 'Cr Condensation Curve')
plt.legend()
plt.tight_layout()
plt.show()
[11]:
# once again we can do a quick sanity check to make sure a cloud is present
show(jpi.mixing_ratio(out['spectrum_output']['full_output'], limit=14, height=600, width=600))
[12]:
# once again we can do a quick sanity check to make sure a cloud is present
show(vpi.all_optics_1d(out['virga_output'], wave_range=[1,2]))
Cloudy vs Clear Spectra
[13]:
opa_mon = jdi.opannection(wave_range=[0.3,15])
hi_res = jdi.inputs(calculation="browndwarf") # start a calculation
teff= 400 # Effective Temperature of your Brown Dwarf in K
grav = 1000 # Gravity of your brown dwarf in m/s/s
hi_res.gravity(gravity=grav, gravity_unit=u.Unit('m/(s**2)')) # input gravity
hi_res.atmosphere(df=out['ptchem_df'])
# now let's add the cloud optical properties
hi_res.clouds(df=out['cld_df'])
df_spec = hi_res.spectrum(opa_mon, calculation='thermal',full_output=True)
wno, fp = df_spec['wavenumber'], df_spec['thermal'] #erg/cm2/s/cm
wno,fp = jdi.mean_regrid(wno,fp, R =200)
[14]:
jpi.show(jpi.spectrum(wno,fp,x_axis_type='log',y_axis_type='log'))
This next spectrum is from the Sonora Elf Owl grid of models which are cloud-free to use as reference. To download the rest of the Sonora Elf Owl grid you can find them separated by L, T, and Y dwarfs.
[15]:
ds_elfowl = xr.load_dataset("spectra_logzz_9.0_teff_400.0_grav_1000.0_mh_0.0_co_1.0.nc")
Now let’s regrid the spectra to R=200 to make sure the two spectra are on the same grid
[16]:
wno_elfowl, fp_elfowl = jdi.mean_regrid(1e4/ds_elfowl['wavelength'].values,ds_elfowl['flux'].values,R=200)
[17]:
fig = plt.figure(figsize=(12,6))
plt.loglog(1e4/wno_elfowl,fp_elfowl, 'k', label = 'Sonora Elf Owl')
plt.loglog(1e4/wno,fp, label = 'Our Cloudy Run')
plt.xlabel('Wavelength [micron]')
plt.ylabel('F$_\\nu$ [erg/cm$^2$/s/Hz]')
plt.xlim(0.5,12)
plt.legend()
plt.show()
For these cloud species you can see their features in the Y and J bands for example in comparison to the cloud-free Sonora Elf Owl model