Climate/03_one_layer_atmosphere_CO2_H2O_model.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<b> One-Layer Atmosphere Model - 2XCO2, Water-Wapor feedback </b><br>\n",
    "Reference: Walter A. Robinson, Modeling Dynamic Climate Systems"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "plt.style.use('seaborn-dark')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Step size\n",
    "dt = 0.01\n",
    "\n",
    "# Set up simulation time\n",
    "tmin = 0\n",
    "tmax = 10\n",
    "t = np.arange(tmin, tmax + dt, dt)\n",
    "n = len(t)\n",
    "\n",
    "# Seconds per year\n",
    "seconds_per_year = 365 * 24 * 60 * 60"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Albedo\n",
    "albedo = 0.30\n",
    "\n",
    "# Emissivity\n",
    "e =  0.85\n",
    "e_0 = 0.75\n",
    "t_0 = 273.15 # K\n",
    "t_1 = 5423 # K\n",
    "default_e = e\n",
    "\n",
    "# Absorption\n",
    "atmospheric_abs = 0.1\n",
    "\n",
    "# Stefan-Boltzmann constant\n",
    "sigma = 5.6696e-8 # W/m^2*K^4\n",
    "\n",
    "# Solar constant \n",
    "solar_const = 1367 # W/m^2\n",
    "default_solar_const = solar_const # W/m^2\n",
    "\n",
    "# Density of water\n",
    "water_density = 1000 # kg/m^3\n",
    "\n",
    "# Depth of the mixed layer\n",
    "depth_mixed_layer = 50 # m\n",
    "\n",
    "# Specific heat capacity of water\n",
    "spec_heat_water = 4218 # J/kg*K\n",
    "\n",
    "# Specific heat capacity of atmosphere\n",
    "spec_heat_atm = 1004 # J/kg*K\n",
    "\n",
    "# Gravity accelleration\n",
    "g = 9.81 # m/s^2\n",
    "\n",
    "# Atmospheric pressure\n",
    "atm_press = 1e5; # Pa (kg/m*s^2)\n",
    "\n",
    "# Atmospheric mass\n",
    "mass_atm = atm_press/g;  # kg/m^2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Heat capacity water\n",
    "heat_capacity_water = water_density * depth_mixed_layer * spec_heat_water # J/K\n",
    "\n",
    "# Heat capacity atmosphere\n",
    "heat_capacity_atm = mass_atm * spec_heat_atm # J/K"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 1: No 2xCO2, 2: 2xCO2 - No Feedbacks, 3: 2xCO2 + Water-vapor feedback\n",
    "models = {\"model_1\": {\"desc\": \"No 2xCO2\", \n",
    "                      \"stronger_sun\": 0, \n",
    "                      \"h2o_feedback\": 0, \n",
    "                      \"T_final\": 0},\n",
    "          \"model_2\": {\"desc\": \"2xCO2 - No Feedbacks\", \n",
    "                      \"stronger_sun\": 1, \n",
    "                      \"h2o_feedback\": 0 , \n",
    "                      \"T_final\": 0},\n",
    "          \"model_3\": {\"desc\": \"2xCO2 + Water-vapor feedback\", \n",
    "                      \"stronger_sun\": 1, \n",
    "                      \"h2o_feedback\": 1 , \n",
    "                      \"T_final\": 0}}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Store the final surface temperature of each model\n",
    "for model_name, model in models.items():\n",
    "    # Initial temperature\n",
    "    ts = 273.15 # K\n",
    "    ta = 273.15 # K\n",
    "    \n",
    "    if model[\"stronger_sun\"]:\n",
    "        # Stronger sun emulate the doubling of CO2 (2xCO2)\n",
    "        solar_const = default_solar_const * 1.02\n",
    "    else:\n",
    "        solar_const = default_solar_const\n",
    "\n",
    "    if model[\"h2o_feedback\"]:\n",
    "        e = default_e\n",
    "\n",
    "    for _ in range(n):\n",
    "        if model[\"h2o_feedback\"]:\n",
    "            e = max(e_0 + 0.1 * t_1 * (1/t_0 - 1/ts), 0)\n",
    "          \n",
    "        # Inflows (Solar to earth, IR down) - Surface\n",
    "        solar_to_earth = (solar_const/4 * (1 - albedo) * seconds_per_year) * (1 - atmospheric_abs)\n",
    "        ir_down = sigma * e * ta ** 4 * seconds_per_year\n",
    "        \n",
    "        # Outflows (IR) - Surface\n",
    "        ir = sigma * ts ** 4 * seconds_per_year\n",
    "        \n",
    "        # Calculate the temperature - Surface\n",
    "        energy_ts = ts * heat_capacity_water + (solar_to_earth + ir_down - ir) * dt\n",
    "        ts = energy_ts/heat_capacity_water\n",
    "        \n",
    "        # Inflows (solar to atmosphere, IR) - Atmosphere\n",
    "        solar_to_atm = (solar_const/4 * (1 - albedo) * seconds_per_year) * atmospheric_abs\n",
    "        \n",
    "        # Outflows (IR down, IR to space) - Atmosphere\n",
    "        ir_to_space = ir * (1 - e) + ir_down\n",
    "        \n",
    "        # Calculate the temperature - Atmosphere\n",
    "        energy_ta = ta * heat_capacity_atm + (ir + solar_to_atm - ir_down - ir_to_space) * dt\n",
    "        ta = energy_ta/heat_capacity_atm\n",
    "    models[model_name][\"T_final\"] = ts"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "No 2xCO2, T=289K\n"
     ]
    }
   ],
   "source": [
    "# Model: No 2xCO2\n",
    "model_1 = models[\"model_1\"][\"T_final\"]\n",
    "print(f\"No 2xCO2, T={round(model_1)}K\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2xCO2 + No feedbacks, T=290K\n"
     ]
    }
   ],
   "source": [
    "# Model: 2xCO2, No feedbacks\n",
    "model_2 = models[\"model_2\"][\"T_final\"]\n",
    "print(f\"2xCO2 + No feedbacks, T={round(model_2)}K\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2xCO2 + Water-vapor feedback, T=291K\n"
     ]
    }
   ],
   "source": [
    "# Model: 2xCO2 + No feedbacks\n",
    "model_3 = models[\"model_3\"][\"T_final\"]\n",
    "print(f\"2xCO2 + Water-vapor feedback, T={round(model_3)}K\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2xCO2 + No feedbacks, ∆T=1.43\n",
      "2xCO2 + Water-vapor feedback, ∆T=2.31\n"
     ]
    }
   ],
   "source": [
    "# Temperature changes in response to 2xCO2\n",
    "delta_t2 = model_2 - model_1\n",
    "delta_t3 = model_3 - model_1\n",
    "print(f\"2xCO2 + No feedbacks, ∆T={round(delta_t2, 2)}\")\n",
    "print(f\"2xCO2 + Water-vapor feedback, ∆T={round(delta_t3, 2)}\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Water-vapor feedback, ∆T=0.89\n"
     ]
    }
   ],
   "source": [
    "# Temperature change in response to the \"Water-wapor feedback\n",
    "delta_H2O = delta_t3 - delta_t2\n",
    "print(f\"Water-vapor feedback, ∆T={round(delta_H2O, 2)}\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Water-vapor feedback, g=0.38\n"
     ]
    }
   ],
   "source": [
    "# Fraction of the total temperature change that is due to the water-wapor feedback (Gain)\n",
    "gain_H2O = delta_H2O/delta_t3\n",
    "print(f\"Water-vapor feedback, g={round(gain_H2O, 2)}\")"
   ]
  }
 ],
 "metadata": {
  "interpreter": {
   "hash": "7b11ec4e63c99d8162a346e44061eab97ff71e663e3c205d951a4ad2f12c3cfa"
  },
  "kernelspec": {
   "display_name": "Python 3.7.10 64-bit (system)",
   "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.7.10"
  },
  "orig_nbformat": 4
 },
 "nbformat": 4,
 "nbformat_minor": 2
}