{ "cells": [ { "cell_type": "markdown", "id": "2d72d361", "metadata": {}, "source": [ "# One-Dimensional Climate Models: Brown Dwarfs w/ Disequilibrium Chemistry at Solar M/H and C/O\n", "\n", "In this tutorial you will learn how to run 1d climate models with the effects of disequilibrium chemistry as was done for the Elf-OWL Grid [Mukherjee et al. 2024](https://ui.adsabs.harvard.edu/abs/2024arXiv240200756M/abstract) (note this should also be cited if using this code/tutorial). \n", "\n", "What you should already be familiar with: \n", "\n", "- [basics of running/analyzing thermal spectra](https://natashabatalha.github.io/picaso/tutorials.html#basics-of-thermal-emission)\n", "- [how to analyze thermal emission spectra](https://natashabatalha.github.io/picaso/notebooks/workshops/ERS2021/ThermalEmissionTutorial.html)\n", "- [how to run a basic 1d brown dwarf tutorial](https://natashabatalha.github.io/picaso/notebooks/climate/12a_BrownDwarf.html)\n", "\n", "\n", "What should have already downloaded: \n", "\n", "1. [Download](https://zenodo.org/record/5590989#.Yzy2YOzMI8a) 1460 PT, 196 wno Correlated-K Tables from Roxana Lupu to be used by the climate code for opacity \n", "2. [Download](https://zenodo.org/record/5063476/files/structures_m%2B0.0.tar.gz?download=1) the sonora bobcat cloud free `structures_` file so that you can have a simple starting guess \n", "\n", "**NEW:**\n", "\n", "3. [Download the .npy](https://doi.org/10.5281/zenodo.10895826) and place them in picaso_refdata folder/climate_INPUTS/661/\n", "\n", "> **_NOTE:_** Tip for getting data from zenodo: pip install zenodo_get then it you can simply retrieve a zenodo posting via the command zenodo_get 10.5281/zenodo.10895826 \n" ] }, { "cell_type": "markdown", "id": "ebb85a88-a88d-497e-a224-a8df6e738a43", "metadata": {}, "source": [ "### First, check that you have downloaded and placed the correlated-k files in the correct folder" ] }, { "cell_type": "code", "execution_count": null, "id": "6d1b911a-9ae9-4284-8d1f-e25227536ebb", "metadata": {}, "outputs": [], "source": [ "import os;import glob\n", "#\n", "glob.glob(\n", " os.path.join(os.environ['picaso_refdata'],'climate_INPUTS','661','*npy')\n", ")\n", "#should see a list of files e.g., \"/data/reference_data/picaso/reference/climate_INPUTS/661/AlH_1460.npy\"" ] }, { "cell_type": "code", "execution_count": null, "id": "5e0676d7-2757-47cb-8e75-047ebb27f80f", "metadata": {}, "outputs": [], "source": [ "import warnings\n", "warnings.filterwarnings('ignore')\n", "import picaso.justdoit as jdi\n", "import picaso.justplotit as jpi\n", "import astropy.units as u\n", "import numpy as np\n", "import matplotlib.pyplot as plt\n", "#%matplotlib inline\n", "from astropy import constants as const\n", "from astropy import units as u\n", "import sys\n", "import pandas as pd\n" ] }, { "cell_type": "markdown", "id": "5de858d5-3e25-4634-a60b-ad339fd9ee06", "metadata": {}, "source": [ "## Setting up Initial Run (highlighting main differences for disequilibrium)" ] }, { "cell_type": "code", "execution_count": null, "id": "3d5079c4-7677-4a03-aa0d-a2850c5cabcc", "metadata": {}, "outputs": [], "source": [ "mh = '+000' #log metallicity\n", "CtoO = '100'# CtoO ratio relative to solar\n", "\n", "ck_db = f'/data/kcoeff_2020_v3/sonora_2020_feh{mh}_co_{CtoO}.data.196'\n", "sonora_profile_db = '/data/sonora_bobcat/structure/structures_m+0.0' #recommended download #2 above\n", "\n", "opacity_ck = jdi.opannection(ck_db=ck_db) # grab your opacities\n" ] }, { "cell_type": "code", "execution_count": null, "id": "8cf8bc9e-b19b-4180-a77a-376b854c7afe", "metadata": {}, "outputs": [], "source": [ "cl_run = jdi.inputs(calculation=\"browndwarf\", climate = True) # start a calculation\n", "\n", "\n", "tint= 700 \n", "grav = 316 # Gravity of your Planet in m/s/s\n", "\n", "cl_run.gravity(gravity=grav, gravity_unit=u.Unit('m/(s**2)')) # input gravity\n", "cl_run.effective_temp(tint) # input effective temperature\n", "\n", "nlevel = 91 \n" ] }, { "cell_type": "markdown", "id": "6f4ad063-2b66-429e-8221-fb34a4077335", "metadata": {}, "source": [ "We recommend starting with Sonora-Bobcat models as an initial guess. " ] }, { "cell_type": "code", "execution_count": null, "id": "8325ae68-7c4e-489e-a99d-3a10cbd39cba", "metadata": {}, "outputs": [], "source": [ "pressure,temp_guess = np.loadtxt(jdi.os.path.join(\n", " sonora_profile_db,f\"t{tint}g{grav}nc_m0.0.dat\"),\n", " usecols=[1,2],unpack=True, skiprows = 1)\n", "\n", "\n", "nofczns = 1 # number of convective zones initially. Let's not play with this for now.\n", "\n", "nstr_upper = 79 # top most level of guessed convective zone\n", "nstr_deep = nlevel -2 # this is always the case. Dont change this\n", "nstr = np.array([0,nstr_upper,89,0,0,0]) # initial guess of convective zones\n", "\n", "# Here are some other parameters needed for the code.\n", "rfacv = 0.0 #we are focused on a brown dwarf so let's keep this as is\n", "print(mh,CtoO,tint)\n" ] }, { "cell_type": "markdown", "id": "b3117928-c7b4-4101-a932-be83c78ccc6a", "metadata": {}, "source": [ "### Setting K$_{zz}$\n", "\n", "We will add one more concept which is the addition of K$_{zz}$ [cm$^2$/s]. K$_{zz}$ is the eddy diffusion constant, which sets the strength of vertical mixing. In `PICASO` we have two options for K$_{zz}$: \n", " \n", " 1. Constant value: sets a constant at every atmospheric layer\n", " 2. Self consistent (see Eqn. 27 and 28 in [Mukherjee et al 2022](https://arxiv.org/pdf/2208.07836.pdf))\n", "\n", "\n", "**New code parameters**: \n", "\n", "0. `diseq_chem=True` : Turns on disequilibrium chemistry\n", "1. `self_consistent_kzz` : (True/False) This solves self consistently for \n", "2. `save_all_kzz` : (True/False) Similar to `save_all_profiles` this saves your intermediate k_zz values if you are trying to solve for a `self_consistent_kzz=True`.\n", "3. `kz` : constant value if `self_consistent_kzz=False`\n", "4. `gases_fly` : **Important**: determines what gases to include in your climate calculation. if you remove one, it will remove the opacity contirubtion from that gas in your climate calculation\n", "5. `chemeq_first` : Converges a chemical equilibrium model first (helpful for convergence)\n", "\n", "**Which of those 6 do I need change change**\n", "\n", "Likely you will only be changing `kz` and/or, for example, playing around with a `self_consistent_kzz` vs a `constant profile`. Unless you are certain, we recommend the following set of `gases_fly` to remain unchanged. \n" ] }, { "cell_type": "code", "execution_count": null, "id": "dd563701-7b34-4b99-8b4d-5f91bdb11970", "metadata": {}, "outputs": [], "source": [ "#following elf-owl lets use a constant value for all pressures\n", "kzval = pressure*0+1e2" ] }, { "cell_type": "code", "execution_count": null, "id": "05142f70-6d40-432a-ab47-0398d272a035", "metadata": {}, "outputs": [], "source": [ "cl_run.inputs_climate(temp_guess= temp_guess, pressure= pressure,\n", " nstr = nstr, nofczns = nofczns , rfacv = rfacv, mh =mh, CtoO = CtoO)\n", "\n", "\n", "gases_fly = ['CO','CH4','H2O','NH3','CO2','N2','HCN','H2','PH3','C2H2','Na','K','TiO','VO','FeH']\n", "\n", "out = cl_run.climate(opacity_ck, save_all_profiles = True, as_dict=True,with_spec=True,\n", " save_all_kzz = False, diseq_chem = True, self_consistent_kzz =False, kz = kzval,\n", " on_fly=True,gases_fly=gases_fly, chemeq_first=False)\n" ] }, { "cell_type": "markdown", "id": "05138391-c05c-43b1-a422-e544d3a3cb89", "metadata": {}, "source": [ "## Compare Diseq and Chemeq Climate Profile \n", "\n", "For the case we chose with very low kzz, and solar M/H the disequilibrium profile and bobcat profiles are identical! " ] }, { "cell_type": "code", "execution_count": null, "id": "71729e76-b6dd-4e38-aa48-65204276da8a", "metadata": {}, "outputs": [], "source": [ "plt.ylim(200,1.7e-4)\n", "plt.semilogy(out['temperature'],out['pressure'],\"r\", label='Elf-OWL Style, Disequilibrium')\n", "plt.semilogy(temp_guess,pressure,color=\"k\",linestyle=\"--\", label='Bobcat, Chemical Equilibrium')" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.11.5" }, "varInspector": { "cols": { "lenName": 16, "lenType": 16, "lenVar": 40 }, "kernels_config": { "python": { "delete_cmd_postfix": "", "delete_cmd_prefix": "del ", "library": "var_list.py", "varRefreshCmd": "print(var_dic_list())" }, "r": { "delete_cmd_postfix": ") ", "delete_cmd_prefix": "rm(", "library": "var_list.r", "varRefreshCmd": "cat(var_dic_list()) " } }, "types_to_exclude": [ "module", "function", "builtin_function_or_method", "instance", "_Feature" ], "window_display": false } }, "nbformat": 4, "nbformat_minor": 5 }