sign_issue.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Think (compute) my way through the problem with the density layers\n",
    "\n",
    "I am having a tough time wrapping my head around this.\n",
    "\n",
    "\n",
    "Lets take isopycnal layers $a$ and $b$, with thicknesses $h_a$ and $h_b$ and tracer concentrations $\\phi_a$ and $\\phi_b$. Together they form a thicker layer $c$, e.g.\n",
    "$h_c = h_a + h_b$ and $phi_c = \\frac{\\phi_a h_a + \\phi_b h_b}{h_a + h_b}$\n",
    "\n",
    "Can it be that $\\phi_c$ decreases in time, even if both $\\phi_a$ and $\\phi_b$ increase, purely to the way that $h_a$ and $h_b$ change?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 140,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 141,
   "metadata": {},
   "outputs": [],
   "source": [
    "def phi_c(phia, phib, ha, hb):\n",
    "    if ha+hb <= 0:\n",
    "        return np.nan\n",
    "    else:\n",
    "        return ((phia*ha)+(phib*hb)) / (ha+hb)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 142,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.9\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "<ipython-input-142-66b70e96e3f1>:29: MatplotlibDeprecationWarning: shading='flat' when X and Y have the same dimensions as C is deprecated since 3.3.  Either specify the corners of the quadrilaterals with X and Y, or pass shading='auto', 'nearest' or 'gouraud', or set rcParams['pcolor.shading'].  This will become an error two minor releases later.\n",
      "  plt.pcolormesh(delta_h_a, delta_h_b, delta_phi_c_results.T,cmap='RdBu_r', vmin=-0.2, vmax=0.2)\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "Text(0, 0.5, 'delta h_b')"
      ]
     },
     "execution_count": 142,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYkAAAELCAYAAAAspXpuAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAAAkSklEQVR4nO3deZRlZXnv8e+vqquq5wG6mbpBOgZFropDM6gZJEoEEi9iTK5D4Er0dlxXjObGFcywjIl3IGrG5cAtCRdNTDAqmNZFHBODWYDScpGpRfsiQjPYtjTQQE/V9dw/9i48XZxT5z3v3nWm/n3W2qvrnLOfvd9NN+etd3peRQRmZmbNjPS6AGZm1r9cSZiZWUuuJMzMrCVXEmZm1pIrCTMza8mVhJmZteRKwsysz0g6S9KdkrZKeleTz98g6ZbyuE7SyamxHZfF6yTMzPqHpFHgu8CZwDbgRuB1EXFHwzkvBrZExE5JZwPviYjTUmI75ZaEmVl/ORXYGhF3RcQ+4Erg3MYTIuK6iNhZvrwBWJca26kFVYL71UKNxrLhfDQbEguUFzeqvMBu3w9gdCQvVqO5cXm/89752GM7ImJNVnDpWC2KPUwnnbuDfbcDexremoyIyYbXa4F7G15vA06b45JvAv45M7atofwmXcYCfoWje10Ms5YOGxvNiluVGXfYeN4X6IrM+wEsXTqeFTexPDduIivuJdd+/QdZgQ32MJ38nfO/+cGeiNgwxynNasmm4wKSzqCoJH6m09hUQ1lJmJl1k4DkBlD7r+xtwLENr9cB9z/lntJzgcuAsyPix53EdsKVhJlZRQLGU7vXDrQ940bgBEnrgfuA1wKvP+h+0nHAVcD5EfHdTmI75UrCzKyioiWRP37TKCKmJF0EfBEYBS6PiNslvaX8/FLg3cDhwIdV3HcqIja0iq1SHlcSZmZVqYPupgQRcQ1wzaz3Lm34+c3Am1Njq3AlYWZWUZ0tiX7jSsLMrKKOBq4HjCsJswqSBytnWZgZtyjzm2hR5hqChZlxAAsW5n29ZMct6uXXmYa2JdHzFdcJOUpWSPqcpG9Lul3Shb0op5lZKwLGpKRj0PS0JVHmGfkQDXlGJG2alWfkrcAdEfFKSWuAOyV9olxybmbWc6p54Lqf9Lq76ck8IwCSZvKMNFYSASxTMc9rKfAQMNXtgpqZzWVYu5t6XUmk5Bn5ILCJYtXgMuA/RURakhQzsy4Y5oHrXo9JpOQZeQVwM3AM8Dzgg5KWP+VC0kZJmyVt3pOwpNHMrC4zU2BTjkHT60oiJc/IhcBVUdgKfB84cfaFImKyXHG4YSH5ScnMzHKMKu0YNL2uJJ7MMyJpnCLPyKZZ59wDvAxA0pHAM4G7ulpKM7M5SMV06JRj0PR0TCIxR8l7gSsk3UrRqrs4Inb0rNBmDbq9biH/fnlxY0vHsuIgf91Cbtz4kvyyVjXMYxK9HrhOyVFyP/CL3S6XmVkqDfFiup5XEmZmw8AtCTMza6pYTDectYQrCTOzijradGjAuJIwM6vIA9dmZjYndzeZDbHc3wK7P5U1737jY3kLTMcy03ZD/pTU8SXjWXG5KcbrIMHIkFYSvV5MZ2Y2BIRG046kq7XfQuFESddL2ivpnbM+u1vSrZJulrS56pO5JWFmVpVgpKZBicQtFB4Cfgt4VYvLnFHXomNXEmZmFUkwMl5bzri2WyhExHZgu6Rfquumrbi7ycysKomR0bQjQbMtFNZ2UJoAviTpW5I2dhDXlFsSZmY10Ejy79yrZ40VTEbEZOOlmsTM3kJhLi+JiPslHQF8WdJ3IuLaDuIP4krCzKwidTYmsSMiNszxecoWCi2V+e6IiO2SrqbovsquJNzdZGZWgxpnN6VsodC8DNISSctmfqZIjnpb5iMBbkmYAfnrDxZmpmLodsrv8cyU32MV0m9np/zOLutEVlwtJJT5b2i2lC0UJB0FbAaWA9OS3gGcBKwGrlaxZmMB8PcR8YUq5XElYWZWkQSjY/V1zCRsofAgRTfUbI8CJ9dWEPqgu6ndopHynJeWC0Nul/Rv3S6jmdncxMjoSNIxaHrakkhZNCJpJfBh4KyIuKccsTcz6x8ieTX1oOl1d1PbRSPA64GrIuIeeHIRiZlZ39AQVxK9bvukLBp5BrBK0tfKxSEXdK10ZmaJ3N00P1IWjSwAXgi8DFgEXC/phoj47kEXKlYWbgRYSm3L483M2pJU68B1P+l1JZGyaGQbxeKTx4HHJV1LMXp/UCVRrlicBFijiU5WJ5p1PXV39pTbzHTYuWm0c6exQn6q8Nyyji1ZmBVXC1HbFNh+0+unSlk08k/Az0paIGkxcBqwpcvlNDObU425m/pKT1sSKYtGImKLpC8AtwDTwGURUWkFoZlZrZS+V8Sg6XV3U9tFI+Xr9wPv72a5zMxSaYi7m3peSZiZDYNB7EpK4UrCzKyiYnbTcM6qdCVhZlbVEC+mcyVhZlaVxyTMBsN4Zuruhem7ih1k6YK8+y1dkHe/sSXjWXHjXY6rdM/li7PixpYsyoqrhzrZmW6guJIwM6uo2JnOlYSZmTVV36ZD/caVhJlZVRIjY8P5dTqcT2Vm1k0CjXoKrJmZNaEh7m4azqcyM+smwcjISNKRdLk22zpLOlHS9ZL2SnpnJ7GdckvChkru1NJux40vzkujPb40L24sMy73fgBjSyby4hbnTWVd0MtU4dS3TiJlW2fgIeC3gFdlxHbELQkzs4pUDlynHAme3NY5IvYBM9s6PykitkfEjcD+TmM75ZaEmVlV6mhMYrWkzQ2vJ8tN02Y029b5tMRrV4ltypWEmVlVnaXl2BERG+a+2lOk7rZZJbYpVxJmZjWoccV1yrbO8xHbVM/HJFJH4iWdIumApNd0s3xmZu1IRe6mlCNByrbO8xHbVE9bEqkj8eV5f0qxzamZWd+pa3ZTyrbOko4CNgPLgWlJ7wBOiohHm8VWKU+vu5ueHIkHkDQzEj97utbbgM8Ap3S3eGZmCWpOy9FuW+eIeJCiKykptopeVxJtR+IlrQXOA36BOSoJSRuBjQBLGc7l8YeS3P1blmT+Npeb8nvJgrx/a+NLMtc7ZMZNLM9bszC+LH/twdjyJVlxuesdxpflpRivgyRGnJZjXqSMxP8lcHFEHJBa/49cTiGbBFijiUqj+WZmnRrWtBy9riRSRuI3AFeWFcRq4BxJUxHx2a6U0MysHe9MN2+eHIkH7qMYiX994wkRsX7mZ0lXAJ93BWFm/cU7082LlFH8XpbPzCyF3JKYP+1G8We9/8ZulMnMrCPedMjMzFqS0IhnN5l1zaCk/M5NwZ07lXV8yXhmXG45K0yBXZw5lXV53lTW8cwpt7VxJWFmZs0JPHBtZmZNeY9rMzNrSYIFeV2B/c6VhJlZRfI6CTMza0l44NrMzFqRKwkzM2vN3U1mXTQoKb9zU3BPrMiMW565TqIHaw8mVi3Nu+eyvHuOLsm7Xy3kloSZmbUioTHPbjIzs6aGdzHdcD6VmVk3zcxuSjlSLiedJelOSVslvavJ55L01+Xnt0h6QcNnd0u6VdLNkjZXfTS3JMzMKqsvwZ+kUeBDwJkUG7PdKGlTRNzRcNrZwAnlcRrwEQ7e+vmMiNhRR3nckjAzq8PISNrR3qnA1oi4KyL2AVcC584651zg41G4AVgp6eh6H6jQ80oioVn1hrI5dYuk6ySd3Itympm1pBG0YDzpAFZL2txwbJx1tbXAvQ2vt5XvpZ4TwJckfavJtTvW0+6mxGbV94Gfj4idks4GJjm4WWV9bNFo3pTUbqf8nliRNzMld0pq7tTZhavypofmTiudWJk/rTT3nguWL8+KG1mSF1cL0cnA9Y6I2NDmarNFB+e8JCLul3QE8GVJ34mIa1MLN1uvWxJtm1URcV1E7Cxf3gCs63IZzczmJIRGR5OOBNuAYxterwPuTz0nImb+3A5cTfE9m63XlURKs6rRm4B/bvaBpI0zzbc9HKixiGZmbdQ7u+lG4ARJ6yWNA68FNs06ZxNwQTnL6XTgkYh4QNISScsAJC0BfhG4rcqj9Xp2U0qzqjhROoOikviZZp9HxCRFVxRrNNH0GmZm86O+FdcRMSXpIuCLwChweUTcLukt5eeXAtcA5wBbgSeAC8vwI4GrJUHx/f73EfGFKuXpdSWR0qxC0nOBy4CzI+LHXSqbmVmyOnM3RcQ1FBVB43uXNvwcwFubxN0F1Dq5J7mSkPRqit/iA/j3iLi6hvs/2awC7qNoVr1+1n2PA64Czo+I79ZwTzOzekkwmrePeL9LqiQkfRj4aeAfyrd+U9LLI+IpNVknEptV7wYOBz5cNqGm2swMMDPrMoF6PcQ7P1JbEj8PPLts4iDpY8CtdRQgoVn1ZuDNddzLzGy+xCFeSdwJHAf8oHx9LHDLvJTIhkruuoUVY3lxyxfnNfkXrlyYF5e73mHloqy4bq93mFi5LCsOur/eYWTZyqy4WohDsyUh6XMUYxArgC2Svlm+Pg24bv6LZ2Y2CFSMSwyhdi2JD3SlFGZmg25IU4XPWUlExL+lXETS9RHxonqKZGY2WEIiRnq9omB+1PVUeR26ZmbD4lAck+iAVzib2SHMU2DNzGwOh/oU2HaGc1jfnpSZ8ZvlC/Ly2awYy4vLTcE9sSI3dXdeT+vEqryppROrMqeyZt4vdxorwMiyVZlxK/PiFudP162FK4k5nV/TdczMBo/qS/DXb5KqPkmvlvQ9SY9IelTSLkmPznweEZVS0ZqZDbrQSNIxaFJbEu8DXhkRW+azMGZmg0mH5jqJBj90BWFm1sIhnJbj1eWPmyV9EvgssHfm84i4av6KZmY2KA7dKbCvbPj5CYqt8GYExT4PZmZ2KFYSEXHhXJ+bmRngtBzzR9JZwF9RbDp0WURcMutzlZ+fQ9GaeWNE3NT1gh7ictctHDae99vV0qXjWXG56xYW5a53yE7BnRe38PAVWXFjK1dmxeWudShi8+45mhmnxXn/bWpTYxbYKt+L7WI71dP2kaRR4EPA2cBJwOsknTTrtLOBE8pjI/CRrhbSzKytckwi5Wh3pQrfi4mxHel1J9qpwNaIuCsi9gFXAufOOudc4ONRuAFYKenobhfUzGwuNa6TqPK9mBLbkeTuJkm/BPwHGjK+RsSfVLk5sBa4t+H1NooNjdqdsxZ4YFb5NlLUqCxlOFc+mlkfSx+4Xi1pc8PryYiYbHhd5XsxJbYjSZWEpEuBxcAZwGXAa4BvVrnxzKWbvDc7o2zKOZT/kScB1mjCWWnNrGtCItLHJHZExIY5Pq/yvZj0fdmJ1KrvxRFxAbAzIv4YeBHFPtdVbZt1nXXA/RnnmJn1TgQHptOOBFW+F2v/vkytJHaXfz4h6RhgP7C+yo1LNwInSFovaRx4LbBp1jmbgAtUOB14JCIemH0hM7NeisQjQZXvxZTYjqSOSXxe0krg/cBNFM96WZUbA0TElKSLgC9STNe6PCJul/SW8vNLgWsopnltpZjq5bUbPbAqewpsXlzuVNbFqxdlxS06vLtTUhd1eSrr6IrDs+IqTYHNvCcTS7LCpjPj6hBAWiMh4VoVvhdbxVYpT3KCv4jYC3xG0ucpBq/3VLnxjIi4huKBG9+7tOHnAN5ax73MzOZL8VVV27WyvxebxVaR2t10fUMB9kbEI43vmZkdymZaEinHoGmX4O8oiilViyQ9n5+MnC+nmO1kZmYBBwawAkjRrrvpFcAbKUbI/7zh/V3A789TmczMBk6d3U39pF2Cv48BH5P0KxHxmS6VycxsoAQw3etCzJN23U3/rdnPMyLiz2e/Z2Z2KBrShkTb7qZlXSmFmdmAG8RB6RTtupv+uFsFsfm3dEF+PsfslN9dX++QN1e+2+sdxg8/LCsue71DZlzu/QBiIm/tyXRmXIzn/ZupQ8Twjkkk/Z8v6RmSvirptvL1cyX94fwWzcxscByItGPQpP56+FHg9yjScRARt1As9zYzO+QV6yQi6Rg0qSuuF0fEN3VwlsOpeSiPmdlAGryv/zSplcQOSU+n/O8g6TXM2s/BzOxQdkgOXDd4K8VeDSdKug/4PvDr81YqM7MBM4A9SUmSKomIuAt4uaQlwEhE7JrfYpmZDY4gmB7SDqfkxXSz3ge8mG7Q5Kb7BlgzkbzT7UEWH56X4mvpEXlTWRcfkZfaeslReVM9J45YnRU3uuqIrsaNLM+bclsl/XZM5C2zivG8fzPT471LFU7AgSFdcp26mO6ZwCn8ZPOKVwLXzlehzMwGSXCIdjfNLKaT9CXgBTPdTJLeA3yqyo0lHQZ8EjgeuBv4tYjYOeucY4GPA0dRpEaZjIi/qnJfM7P5MKzdTanrJI4D9jW83kfx5V7Fu4CvRsQJwFfL17NNAb8TEc8CTgfeKumkivc1M6tdseq6/TFoUjua/xb4pqSrKVpW5wEfq3jvc4GXlj9/DPgacHHjCeWerQ+UP++StIVif4s7Kt7bzKw2M4vphlFSSyIi/gfFHqo7gYeBCyPif1W895FlJTBTGcw5IifpeOD5wDdafL5R0mZJm/dwoGLRzMzSRcD+A5F0VCHpMElflvS98s+mMzUknSXpTklbJb2r4f33SLpP0s3lcU67eyZPWYmIm4CbUs8vC/QVivGE2f6gw+ssBT4DvCMiHm1RvkmKtRys0cRwVulm1qeCA91pScx0019Sfvm/i1k9MJJGgQ8BZwLbgBslbYqImR6Yv4iID6TeMG9eY6KIeHmrzyT9UNLREfGApKOB7S3OG6OoID4REVfNU1HNzLJ1sbupbTc9cCqwtVzfhqQry7isbvp5rSTa2AT8Z+CS8s9/mn2CigUZfwNs8ZqMnxhV+3OaOWph/jqJpWsy1zsck5f2efHRefP6Fx+VF5e93uHwo/PiVq3JitPyvHJOL8xbs9CLdRL7lfe1tHt/DzsQOlsnsVrS5obXk2VPSIqDuuklNeumXwvc2/B6G3Baw+uLJF0AbKaYGHTQrNLZellJXAL8o6Q3AfcAvwog6Rjgsog4B3gJcD5wq6Sby7jfj4hrelBeM7OmOmxJ7IiIDa0+rKGbvtmvkTOF+wjw3vL1e4E/A35jrov1rJKIiB8DL2vy/v3AOeXP/07zBzYz6yt1jUnU0E2/DTi24fU64P7y2j9suNZHgc+3K0/+VmVmZgZ0b3YTP+mmhxbd9MCNwAmS1ksap9j7ZxNAWbHMOA+4rd0Ne9ndZGY2FIKubSjUtps+IqYkXQR8ERgFLo+I28v490l6HkV3093Ab7a7oSsJM7MadGNr0pRu+vL1NcBTxm4j4vxO7+lKwsysomFece1KYgAdmZm2+8iFY9n3zJ3KuvToFXlxa/Omei46+sisuAVHrMuKGz282SSUBEvyUppPL8r77zm9cHleXGbaboDd+/NyZ++eyovb28tc3QEHhnRrOlcSZmYVBbDflYSZmTXj7iYzM2stgmm3JMzMrJmgO7ObesGVhJlZDdzdZGZmTRUtCVcSZmbWxExajmHkSmIArV2U99e24ml5c+UBVhy3Mitu2bF56xYWr81Lwb3gqOPy4taszYqbXpKXmjx/vUNe3J7pvDyZu/fk7/KYu95hz1Tel23uuoy6uLvJzMyaiu7tTNd1riTMzKoa4hXXPUsVnrqhd3nuqKT/K6lt7nMzs24Likoi5Rg0vdxPYmZD7xOAr5avW3k7sKUrpTIz61CEK4n5cC7FRt6Uf76q2UmS1gG/BFzWnWKZmXUmCPZNTScdg6aXYxIpG3oD/CXwu8Ccu6pL2ghsBFjKaI3FNDNrY4jHJOa1kqi6obekXwa2R8S3JL10rnMjYhKYBFijiYH42zpsPK8yW7tqUVbcyuPzpk8CLF+fNyV16fpj25/UxIJj1ufFHfm0rLgDSw7PiptenJfye9/oRFbc45nTPJ/YnzeVdff+/P+Vdu2byop7bF9eWXftzbtfHWbGJIbRvFYSNWzo/RLgP0o6B1gILJf0dxHx6/NUZDOzjsUQtyR6OSbRdkPviPi9iFgXEcdTbOb9L64gzKwfdWPgOnVWqKTLJW2XdFtOfKNeVhKXAGdK+h5wZvkaScdIesrerGZm/Wo6gr1T00lHRamzQq8AzqoQ/6SeDVynbujd8P7XgK/Ne8HMzDJ0qbvpXOCl5c8fo/hOvHj2SRFxraTjc+MbecW1mVlFXRyTSJ0VWlu8Kwkzsxp0kLtptaTNDa8ny9mZQPVZoXVzJWFmVlHQ0aD0jojY0PJa1WeFzqXjeFcSPfT0JWNZcYedkDc3f9Uz8tJhA6x4Rt66hbHjnpEVN3JE3nqHqWVrsuIOLM5L+b1rb96c/scy5/TnrpN4eHfumoX8tQePZD5j/n+b/LTmVXWxu2lmVugltJgVWnd8L2c3mZkNhQD2TR1IOipKmhUq6R+A64FnStom6U1zxc/FLQkzs6qiO8n7UmeFRsTrOomfiysJM7OKnJbDzMxaioApVxJmZtaMWxJmZtbaECf4cyVRg1Hlxa1ftzwr7oiT89JvH/bcZ2bFAYz/9HPzAo/Imzo7tbzZWqL2HpvOS7/+yGP78+63L29K6o4n9mXF7dydV87c6ahV0m8/8kReWR/OjMudOluHmU2HhpErCTOzioY5VbgrCTOzGoQrCTMzayYCpl1JmJlZc0GkJ/gbKD1Ly9HBDksrJX1a0nckbZH0om6X1cxsTgEHpqaTjkHTy9xNqTsk/RXwhYg4ETgZ2NKl8pmZJQkgptOOQdPLSuJcip2RKP981ewTJC0Hfg74G4CI2BcRD3epfGZmySIi6Rg0vRyTSNkh6aeAHwH/R9LJwLeAt0fE410sZ1svWLkwK27t6euy4o449TlZcQufk99TN7U6b73DE+Mrs+Ie2p2XLXPn7rz1B9sf35sVtyNzTv9DmeV86LG8uB9nxj2SWU6AXXsyU4Vnxu3v4ToJhnjgel5bEpK+Ium2Jse5iZdYALwA+EhEPB94nBbdUpI2StosafMeepdX3swORUFMpx2DZl5bEjXssLQN2BYR3yhff5oWlUS5/d8kwBpNDN7fhJkNrGJMYji/dno5JjGzQxK02CEpIh4E7pU0k0/iZcAd3SmemVmigAMHppOOQdPLMYlLgH8sd0y6B/hVKHZYAi6LiJkNNN4GfELSOHAXcGEvCmtmNpdhbUn0rJLoYIelm4GWm4abmfVaRHjg2szMWuvGFNgOFiFfLmm7pNtmvf8eSfdJurk8zmkW38hpOWrwnFOPyYo79hUvzoqb2HBmVtzeI/NThT+QmUr7vu278+J27cmK++FjeVNZH3w4737bd+Xd78eZ5dyVOeV27+7uTyvdlxk7tT+v335qf29nNXZpodzMIuRLJL2rfH1xk/OuAD4IfLzJZ38RER9IvaFbEmZmFc0k+Es5Kmq7CLkoT1wLPFT1ZuBKwsysuoDpqemko6KDFiEDzRYht3ORpFvKLqmm3VWNXEmYmVUWTEfaAayeWfhbHhsbr1TDIuS5fAR4OvA84AHgz9oFeEzCzKyiDhfT7YiIljM2a1iE3LqcET9suNZHgc+3i3FLwsysqqBbaTnaLkKeS1mxzDgPuK3VuTNcSZiZ1aBLA9eXAGdK+h5wZvkaScdIumbmJEn/AFwPPFPStnLRMsD7JN0q6RbgDOC3293Q3U1mZhVFBNNdSLnRwSLk17WIP7/Te7qSaHDKqryU38/+L23XozS14IzXZ8X9YHpZVtyddz+SFQfw3R/nZWf/3oOPZcVt2/lEVtzOR/LWO+x5PG/9wZ7MdQu5awj27c6739S+vHUZB/bmrXMBOLAvL3YqM256f35a8zoM64prVxJmZjWI6eHcosCVhJlZVRGuJMzMrLnAlYSZmbUSEAdcSZiZWTMxzfRUbwfO54srCTOzGri7qWaSDgM+CRwP3A38WkTsbHLebwNvplj5fitwYUTkzXNs43UffXNW3IM//5tZcZ/b8qOsuH/dck9W3AP378qKA3h0Z960xCceyZsCu29XXgLLfY/nTfOd2pM3xTd3umbu1NJh/W110A3zmEQvV1zP5EU/Afhq+fogktYCvwVsiIhnA6PAa7taSjOzdqJoSaQcg6aXlURSXnSK1s4iSQuAxcD98180M7NOBNPTB5KOQdPLMYmD8qJLekpe9Ii4T9IHgHuA3cCXIuJLzS5WptvdCLCU0fkrtZnZbF4nkUfSV4Cjmnz0B4nxqyhaHOuBh4FPSfr1iPi72edGxCQwCbBGE8O5Pt7M+lJE9DwtyHyZ10qihrzoLwe+HxE/KmOuAl4MPKWSMDPrpWFtSfRyTCIlL/o9wOmSFksSRfbDLV0qn5lZmrK7yQPX9WqbFz0ivgF8GriJYvrrCGWXkplZ/xjeSqJnA9cd5EX/I+CPOrn20174HC7dvLnjMo0//zc6jgHgv789L87MhkKxfen87yfRC15xbWZWVcTQLnR0JWFmVlXEQK6BSOFKwsysomB4s8D2cuDazGw4dGl2k6TDJH1Z0vfKP1c1OedYSf8qaYuk2yW9vZP42VxJmJlV1rXZTW1z3gFTwO9ExLOA04G3Sjqpg/iDuJIwM6tBlyqJtjnvIuKBiLip/HkXxdqytanxsyli+DJYSPoR8IMu3W41sKNL9+qmYXyuYXwm8HNV9bSIWFPlApK+QFHeFAuBxu0OJsu0Qin3eTgiVja83hkRLbuMJB0PXAs8OyIe7TQehnTguupfeCckbY6IDd26X7cM43MN4zOBn6sfRMRZdV2ras67hussBT4DvCMiHs0tz1BWEmZmg6qGnHdIGqOoID4REVc1fJQU38hjEmZmg6Ntzrsyz93fAFsi4s87jZ/NlUR1w5pLahifaxifCfxch5K2Oe+AlwDnA78g6ebyOGeu+LkM5cC1mZnVwy0JMzNryZWEmZm15EoikaSzJN0paaukp6xSVOGvy89vkfSCXpSzEwnP9IbyWW6RdJ2kk3tRzk61e66G806RdEDSa7pZvlwpzyXppWUf9O2S/q3bZexUwr/BFZI+J+nb5TNd2ItyHtIiwkebAxgF/h/wU8A48G3gpFnnnAP8MyCKpfDf6HW5a3imFwOryp/P7vdnSn2uhvP+BbgGeE2vy13T39dK4A7guPL1Eb0udw3P9PvAn5Y/rwEeAsZ7XfZD6XBLIs2pwNaIuCsi9gFXUixvb3Qu8PEo3ACsLOch96u2zxQR10XEzvLlDcC6LpcxR8rfFcDbKOaRt50n3idSnuv1wFURcQ9ARPT7s6U8UwDLymmdSykqianuFvPQ5koizVrg3obX2/hJLpROzuknnZb3TRQtpX7X9rkkrQXOAy7tYrmqSvn7egawStLXJH1L0gVdK12elGf6IPAs4H6KLYzfHhHDuQVcn/KK6zRq8t7sucMp5/ST5PJKOoOikviZeS1RPVKe6y+BiyPiQPEL6kBIea4FwAsptgVeBFwv6YaI+O58Fy5TyjO9ArgZ+AXg6cCXJX09KqSZsM64kkizDTi24fU6it9sOj2nnySVV9JzgcuAs6PYl7zfpTzXBuDKsoJYDZwjaSoiPtuVEuZJ/Te4IyIeBx6XdC1wMtCvlUTKM10IXBLFoMRWSd8HTgS+2Z0imrub0twInCBpvaRx4LUUy9sbbQIuKGc5nQ48EhEPdLugHWj7TJKOA64Czu/j30Zna/tcEbE+Io6PiOOBTwP/tc8rCEj7N/hPwM9KWiBpMXAaRZrofpXyTPdQtIyQdCTwTOCurpbyEOeWRIKImJJ0EfBFihkZl0fE7ZLeUn5+KcUsmXOArcATFL8B9a3EZ3o3cDjw4fK37qno86ycic81cFKeKyK2lCmrbwGmgcsi4rbelXpuiX9X7wWukHQrRffUxRExjGnR+5bTcpiZWUvubjIzs5ZcSZiZWUuuJMzMrCVXEmZm1pIrCTMza8mVhJmZteRKwvqapPdIemfqOZLeKOmYDu9xxaCkCzfrNlcSNmzeCHRUSZhZa64krO9I+oNyI5qvUKRhmHn/6ZK+UGY4/bqkE2fFvYYiL9Mnyo13Fkl6t6QbJd0maVKtM/r9XLmx0l1ztSokLZX0VUk3SbpVUrM05GZDw5WE9RVJL6TI4fN84NXAKQ0fTwJvi4gXAu8EPtwYGxGfBjYDb4iI50XEbuCDEXFKRDybIjPqL7e49dEUWW5/GbhkjiLuAc6LiBcAZwB/NkfFYzbwnLvJ+s3PAldHxBMAkjaVfy6l2CnvUw3fyRMJ1ztD0u8Ci4HDgNuBzzU577PlPgV3lInkWhHwPyX9HEV+pLXAkcCDCWUxGziuJKwfNUsoNgI8HBHPS72IpIUUrY0NEXGvpPcAC1ucvrcxdI7LvoFiG80XRsR+SXfPcU2zgefuJus31wLnleMJy4BXApSbzHxf0q8ClCnZT24SvwtYVv488+W9o2yJ1DGDaQWwvawgzgCeVsM1zfqWWxLWVyLiJkmfpNiN7AfA1xs+fgPwEUl/CIxR7In87VmXuAK4VNJu4EXARym2vbybYv+Cqj4BfE7S5rKM36nhmmZ9y6nCzcysJXc3mZlZS+5uMmtC0nOAv5319t6IOK0X5THrFXc3mZlZS+5uMjOzllxJmJlZS64kzMysJVcSZmbW0v8HPBRyh2/CG0oAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# both concentrations on each layer increase\n",
    "phi_a_0 = 0.2\n",
    "phi_a_1 = 0.4\n",
    "\n",
    "phi_b_0 = 0.7\n",
    "phi_b_1 = 1.0\n",
    "\n",
    "\n",
    "ha = 0.1\n",
    "hb = 0.9\n",
    "\n",
    "# only works if the thicker layer has the higher concentration? That is similar to my problem\n",
    "\n",
    "# the concentration values in the upper layer cannot approach the lower layer either...that is also very likely in my example...although I have to check that...\n",
    "\n",
    "\n",
    "max_change = np.max([ha,hb]) * 1\n",
    "print(max_change)\n",
    "\n",
    "delta_h_a = np.linspace(-ha,max_change, 20)\n",
    "delta_h_b = np.linspace(-hb,max_change, 30)\n",
    "delta_phi_c_results = np.ones([20,30])\n",
    "\n",
    "for ii,dha in enumerate(delta_h_a):\n",
    "    for jj,dhb in enumerate(delta_h_b):\n",
    "        delta_phi_c = phi_c(phi_a_1, phi_b_1, ha+dha, hb+dhb) - phi_c(phi_a_0, phi_b_0, ha, hb)\n",
    "        delta_phi_c_results[ii, jj] = delta_phi_c\n",
    "        \n",
    "plt.pcolormesh(delta_h_a, delta_h_b, delta_phi_c_results.T,cmap='RdBu_r', vmin=-0.2, vmax=0.2)\n",
    "plt.colorbar()\n",
    "plt.xlabel('delta h_a')\n",
    "plt.ylabel('delta h_b')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Am I convinced? Meh...."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Plug in realistic values"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 143,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "<ipython-input-143-87e2e54775b9>:20: MatplotlibDeprecationWarning: shading='flat' when X and Y have the same dimensions as C is deprecated since 3.3.  Either specify the corners of the quadrilaterals with X and Y, or pass shading='auto', 'nearest' or 'gouraud', or set rcParams['pcolor.shading'].  This will become an error two minor releases later.\n",
      "  plt.pcolormesh(delta_h_a, delta_h_b, delta_phi_c_results.T,cmap='RdBu_r', vmin=-50, vmax=50)\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.lines.Line2D at 0x2af18d63cc10>"
      ]
     },
     "execution_count": 143,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYMAAAELCAYAAAA7h+qnAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAAAfuElEQVR4nO3dfZRdVZ3m8e9TlRfyCkJeSQXJtFEnIKOSjth2+zKiRAWCimulRxtsXZOWhT326nbZxPRydFyZRWs3dqsDTkkraKORHkUiLYOEbsd2NRhiRCFAJBICRULeIK+QSlL1mz/OKXNT3Lr31KlTde+59/msdVbdu8/Z5+xNtH61X87eigjMzKy9dTS6AGZm1ngOBmZm5mBgZmYOBmZmhoOBmZnhYGBmZpQ0GEhaKmmzpC2Srml0eczMyk5le89AUifwa+BtQA9wP/CHEfFwQwtmZlZiZWwZLAG2RMTjEXEUWAMsa3CZzMxKbVyjC5DDPOCpiu89wOsGXyRpBbACYBw6/zTGj03p2sBp4ztz5Zs8a2rN80/2TQHgrM7DJ6WPf8mMXM872nlKrnz7eo/lyrf3QG+ufEcOHMiVL/qO58pnQ4sX9u6JiJkjucd8TYoj9Ge6dg9H74qIpSN5XlHKGAxUJe1FfV0R0Q10A8zUxHgvc0e7XG3jktnTc+Vb/Kdvqnn+qgNLALhh+vqT0mdefmWu520/7ZW58n1/8+5c+b5595Zc+R5ed1eufL0Hn82Vz4Z27IGvbxvpPY7QT9bfN/+bbfn+0hkFZQwGPcD8iu9dwPYGlcXM7CQCOqv9yVpNEw3ZljEY3A8slLQAeBpYDvyXxhbJzCwhYEJHxmjQN6pFGZbSBYOIOC7po8BdQCfwtYjY1OBimZkBAy2DrE2D5lG6YAAQET8EftjocpiZvYiG0U3UREoZDMzMmpVbBmZmNrwB5CbiYGDDNu3M2u8LDGXynDNqnu84Mq7qdX2TX5LreQeOZpvrPdjTz76QK9+hfUdy5fMU0VYjtwzMzNqdgPEOBmZm7U0eQDYzM/AAsplZ2/MAspmZeWqpmZkl3DIwM2tz0jDWJmoiDgZtbOq4fHsbTZk1JVe+8TNm1zyvpydUva4/53sGu3cfzZVv297D9S+q4uCz+3Lls9biMQMzM0N+6czMzKCcLYMy7oFsZta0kpfOlOnIdj91SvqFpDvS76dLulvSY+nPfP2ogzgYmJkVaGBzmyxHRh8DHqn4fg1wT0QsBO5Jv4+Yg4GZWYEGBpCzHHXvJXUB7wJurEheBtycfr4ZuKyIcnvMwMysYMMYQJ4haUPF9+6I6K74/nfAJ4BpFWmzI2IHQETskDRrJGUd4GDQxs48Jd8//9S50+pfVEXnGXNrntf4zqrXHRk3Kdfzdh5+Lle+HTmXsH5+79O58llrkaAjezDYExGLq99HFwO7IuLnkt5cUPGG5GBgZlYooWKmE70BuFTSO4FTgOmS/hHYKWlu2iqYC+wq4mEeMzAzK5Kgo1OZjloiYmVEdEXE2cBy4F8i4gPAWuDK9LIrgduLKLZbBmZmBZKgY0LnaD7iWuBWSR8GngTeV8RNHQzMzIqk+n/1D1dE/Bj4cfp5L/DWQh+Ag4GZWeHUUb4eeAcDM7MCKR0zKBsHAzOzghU0m2hMORi0sXmT8v3zT5l7Rq58nWfMqXle4/an1516UvqB3v5cz9t1uDdXvoPP5XvP4PCup3LlsxYjoU53E5mZtTUJOsc7GJiZtTnR4ZaBmVmbk8cMzMzanhwMzMwMcDeRmVm7k+QBZCuX08/ItzT0lDn5ppb2T669O190Hq563f7evlzPe3LP87nyHdp3JFe+vqP5pqRaixGeWmpmZn4D2czMVNh+BmOqYW0ZSe+TtElSv6TFg86tlLRF0mZJF1Wkny/pwfTcF6Xs2wmZmY0Fpd1EWY5m0sjSPAS8B/hJZaKkRSQbOZwDLAWulzSwOPgNwApgYXosHbPSmpllVMTmNmOtYd1EEfEIJCPvgywD1kREL7BV0hZgiaQngOkRcW+a7xvAZcCdY1VmM7N6ktlEo7q5zahoxjGDecB9Fd970rRj6efB6VVJWkHSimAq5fuHMbOS8ktnLyZpHVBtqcpVETHUvp3V/itGjfSqIqIb6AaYqYlDXmdmVihPLX2xiLgwR7YeYH7F9y5ge5reVSXdcpo6d2qufBNnzciVr957BnTsrHrd3meP5Xretr2Hc+U7tGdPrnxmCZVyp7NmLPFaYLmkiZIWkAwUr4+IHcBBSReks4iuAIZqXZiZNUSy01lHpqOZNGzMQNK7gS8BM4F/lvRARFwUEZsk3Qo8DBwHro6IgVdQrwJuAiaRDBx78NjMmow3txmWiLgNuG2Ic6uB1VXSNwDnjnLRzMzyk+gY34xzc2orX4nNzJqZQJ3lm8HoYGBmViC5m8jMzBB0lHA2kYNByc2cmL85Oi3n1NJxs7rqX1TF0QnTap7v7+iset2OQ8/let7O5/ItKX1495O58pkNcMvAzKzNyQPIZmaWLGHtloGZWXvzchRmZgY03dvFWTgYmJkVSCrn2kQOBmZmBXM3kZlZu/NsImuEM08ZnzvvlLmn58rXeUa1LSrqO3i0v+b5vv6oet2Og0fyPe/ZvO8ZPJUrnxmkU0tLuBxF+doyZmZNbvDG90Mdde8jzZf0r5IekbRJ0sfS9NMl3S3psfRnnc1C6nMwMDMrkooLBiTL+P9FRPxH4ALgakmLgGuAeyJiIXBP+n1EHAzMzAqVzCbKctQTETsiYmP6+SDwCMne78uAm9PLbgYuG2mpPWZgZlYgjdJLZ5LOBl4D/AyYne7+SETskDRrpPd3MDAzK9LwZhPNkLSh4nt3RHS/+JaaCnwX+LOIOJDs/FssBwMzsyJJqCPzbKI9EbG49u00niQQ3BIR30uTd0qam7YK5gK78hc44WBQcvMm5f8nnDIn39TSmD4zV779vX01zx+P6tdt2/N8rucdyDm11GzEsgeDmpQ0Af4BeCQirqs4tRa4Erg2/Xn7SJ/lYGBmVihBcctRvAH4I+BBSQ+kaZ8kCQK3Svow8CTwvpE+yMHAzKxIBe6BHBE/Te5Y1VsLeUjKwcDMrEgSjJvQ6FIMm4OBmVmBhFctNTMzUdgA8lhyMDAzK5QcDMzMDHcT2dibNmdK7rxT5uV7X6B/cr4FEnfvO1bz/LG+ZOnq3YdPvm7b3sO5nndo19O58pmNiNwyMDMzCY33bCIzszZX6EtnY8bBwMysSJ5NZGZmyX4GDgZmZuZuIjOzNqcO5OUobKxN75qeO2/nGXNz5eubnG/p6x09z9Y8f7QvWcN6x6Hek9J37c23hPXh3U/lymc2IqKULYOGlVjS5yU9KulXkm6TdFrFuZWStkjaLOmiivTzJT2YnvuiRmO7HzOzERBCnZ2ZjmbSyPB1N3BuRJwH/BpYCSBpEbAcOAdYClwvaeC/2g3ACmBheiwd60KbmdU0MJsoy9FEGtZNFBE/qvh6H3B5+nkZsCYieoGtkrYASyQ9AUyPiHsBJH0DuAy4s96zDp8+i59d9LECS988njpjcu68kzfm6+7hNz/PlW3vC8drnt/+XNId9JV1j52Uvmv/kVzPm7bovFz5phyv/aa0ta5nHvh6AXcp5xvIzdKx9SFO/FKfB1R29vakafPSz4PTq5K0QtIGSRv6Cy6smVkt6ujIdDSTzC0DSe8Bfh8I4KcRcVuGPOuAOVVOrYqI29NrVgHHgVsGslW5PmqkVxUR3UA3wExNjNd9++/rFbeU/uvFC3PnPeeS5bnyxRvfkSvfD35dewB5oEXwkQtPrtOX7tqc63k7fvqrXPmO7N+dK58ZkKxN1Dm+0aUYtkzBQNL1wMuAb6dJfyLpwoi4ula+iLiwzn2vBC4G3hoRA7/Ye4D5FZd1AdvT9K4q6WZmTUSg5vqrP4usLYM3kQz2BoCkm4EHR/JgSUuBvwTeFBGVcwfXAt+SdB1wJslA8fqI6JN0UNIFwM+AK4AvjaQMZmajIVo4GGwGzgK2pd/nA/na4Cd8GZgI3J3OEL0vIj4SEZsk3Qo8TNJ9dHVE9KV5rgJuAiaRjDHUHTwui/mT8jUrp86dlvuZ42YOOeRS057evvoXVbHjYO2B4KPpEtaDrzv43Au5nufuHmsI0XotA0k/IOmXPxV4RNL69PvrgH8fyYMj4mU1zq0GVldJ3wCcO5LnmpmNLiXjBiVTr2XwN2NSCjOzVtJkM4WyqBkMIuL/ZbmJpHsj4vXFFMnMrLxCIjrKt9JPUSU+paD7mJmVX6uNGQzDkPP9zczaS2tPLTUzs4xaeWppPeUbOm8y8ybl+6eYOm9m7mfqJdVeDq9vf86ppdv21F6K+sixvqrX7a+Tz6zptHEw+KOC7mNmVm5q4YXqJL1H0mOS9ks6kL4JfGDgfEQ8NHpFNDMrl1BHpqOZZG0ZfA64JCIeGc3CmJmVn1rvPYMKOx0IzMwyaNHlKN6Tftwg6TvA94HfblAbEd8bvaKZmZVRa04tvaTi8/PA2yu+B+BgYGY2WKsFg4j447EqiJlZS2jz5ShshOaclm9Fj6ld+d8z6J+Sbw/kXQfz7RH8+O5DNc+/cLSv6nWHdz+Z63lmDVPgqqXp3i9/D3QCN0bEtYXdvEL52jJmZk0tHTPIctS7k9QJ/C/gHcAi4A8lLRqNUjsYmJkVrMD3DJYAWyLi8Yg4CqwBlo1GmTN3E0l6F3AOFSuURsT/GI1CmZmVWvYB5BmSNlR8746I7orv84CnKr73kGwuVrhMwUDSV4DJwFuAG4HLgfWjUSAzszILicg+ZrAnIhbXOF/tRqOySnTW8PV7EXEF8FxEfAZ4Pck+yGZmVimCvv5sRwY9nPy7tgvYPhrFzhoMBnYkf17SmcAxYMFoFMjMrOwi45HB/cBCSQskTQCWA2uLL3H2MYM7JJ0GfB7YSFKPG0ejQO1qete0XPkmzJ6b+5l9U/NNS92x47lc+fbsfaHm+WPH+6ted+iZJ3I9z6wRAsj2R3+Ge0Ucl/RR4C6SqaVfi4hNxdz9ZJkXqouIXuC7ku4gGUQ+MhoFMjMru4jiuvUj4ofADwu74RCydhPdO/AhInojYn9lmpmZJQZaBlmOZlJvobo5JFObJkl6DSdGtqeTzC4yM7NKAX1N9os+i3rdRBcBHyQZwb6uIv0g8MlRKpOZWakV2U00VuotVHczcLOk90bEd8eoTGZmpRVAf6MLkUO9bqI/r/Z5QERcNzjNzKzdlbBhULebKN98RzOzNtZsg8NZ1Osm+sxYFaRVTOrMt3Rt3vcMxs2clysfwHNH8zVme/bXfl9gKPv3Pl/zfJ+qX3f8SO2lr82aSUQ5xwwyTS2V9HJJ90h6KP1+nqS/Gt2imZmVU19kO5pJ1vcMvgqsJFmGgoj4Fclr0WZmViF5zyAyHc0k6xvIkyNivU5eie/4KJTHzKz0muvXfDZZg8EeSb9DWkdJlwM7Rq1UZmYl1nIDyBWuBrqBV0p6GtgKfGDUSmVmVmJN1gOUSaZgEBGPAxdKmgJ0RMTB0S2WmVk5BUF/CTuKMr90NigdGNlLZ5I+S7KXZz+wC/hgRGxPz60EPgz0Af8tIu5K088HbgImkazi97Fosjlc8yeNz5VvateMXPk6Z3blygewvzff1NJte2pPER3KgT37a19w2pTkun2Hc93frCkE9JXwFeR6s4mmpcdi4CqSRevmAR8BFo3w2Z+PiPMi4tXAHcCnACQtIpmpdA6wFLheUmea5wZgBbAwPZaOsAxmZoUKBt41qH80k0wvnUn6EfDage4hSZ8G/mkkD46IAxVfp3BiAH4ZsCbdP2GrpC3AEklPANMj4t60DN8ALgPuHEk5zMyK1nLdRBXOAo5WfD8KnD3Sh0taDVwB7AfekibPA+6ruKwnTTuWfh6cPtS9V5C0IphK51CXmZkVrtn+6s8iazD4JrBe0m0kf8G/G7i5XiZJ64A5VU6tiojbI2IVsCodI/go8N85sWdCpaiRXlVEdJPMgGKmJpbwn8bMymjgpbOyyTqbaLWkO4E/SJP+OCJ+kSHfhRnL8S3gn0mCQQ8wv+JcF7A9Te+qkm5m1jQi4FizrTWRQdaWARGxEdhY1IMlLYyIx9KvlwKPpp/XAt+SdB1wJslA8fqI6JN0UNIFwM9Iupe+VFR5zMyKEfS1astglFwr6RUkU0u3kcxQIiI2SboVeJhkyYurI6IvzXMVJ6aW3okHj82sybR0N9FoiIj31ji3GlhdJX0DcO5olmukzpqc8z2DefneM+ibckaufAC7nz9a/6IqHt+db0npQ89srXl+0uSXAfBCnevMmlpJ3zNoZMvAzKzluGVgZmYAHjMwM2t3LT+byMzM6guab+OaLBwMzMwKVsKGgYOBmVmRPIBsAJw+e0qufJPnzs6Vr39qvimpAD07821LsSfnEtaHdj5R8/yErq5M15k1tYC+Em515mBgZlagAI45GJiZtTd3E5mZGUTQ75aBmVl7CzybyMzMKGc3Ub09kM3MbBiSlkFkOkZC0uclPSrpV5Juk3RaxbmVkrZI2izpoiz3czAwMyvQwHIUWY4Ruhs4NyLOA34NrASQtAhYDpwDLAWul1R37193ExVsete0XPnGzTkrV779ffn3d+7Z/0KufAdyvmcQ/X31rsh4nVlzG4tuooj4UcXX+4DL08/LgDUR0QtslbQFWALcW+t+DgZmZgWKxux09iHgO+nneSTBYUBPmlaTg4GZWZGG9wbyDEkbKr53R0T3wBdJ64A5VfKtiojb02tWkewKectAtuqlqs3BwMysQMGwgsGeiFg85L0iLqyVWdKVwMXAWyN+2xzpAeZXXNYFbK9XEA8gm5kVKNKWQZZjJCQtBf4SuDQiKgfy1gLLJU2UtABYCKyvdz+3DMzMChQER4+PySbIXwYmAndLArgvIj4SEZsk3Qo8TNJ9dHVE1J2V4WBgZlakMVq1NCJeVuPcamD1cO7nYDCE0yfkm7J56ktfkivfuNnz619Uxf7e/NMwH991ON8zd+7K/UyzVjfMMYOm4WBgZlag8H4GZmYGDgZmZm2vP4LesRlALpSDgZlZwdwyMDNrcx4zMDMzgEasTTRiDgZmZgUKRv52cSM4GAxh/qTxufJNPWt2rnz902blyrfr8NFc+QAe330oV75DO7fmfqZZq3M3kZmZEcDR4+Xbk8PBwMysSOFuIjOztuflKMzMjAg4XsJg0PD9DCR9XFJImlGRtlLSFkmbJV1UkX6+pAfTc19Uum6rmVmzGGgZjPZ+BkVraDCQNB94G/BkRdoiYDlwDrAUuF7SwBKiNwArSDZrWJieNzNrHmO0uU3RGt1N9AXgE8DtFWnLgDUR0QtslbQFWCLpCWB6RNwLIOkbwGXAnaNRsLOmTciVb+q8mbny9U3Nl69n9/P1LxrC3pxLWD+/t+4OemZtaww3tylUw4KBpEuBpyPil4N6e+YB91V870nTjqWfB6cPdf8VJK0IppJvbwIzs+HyewZVSFoHzKlyahXwSeDt1bJVSYsa6VVFRDfQDTBTE8v3L2NmpRUOBieLiAurpUt6FbAAGGgVdAEbJS0h+Yu/ctuvLmB7mt5VJd3MrGlEQH8Jg0FDBpAj4sGImBURZ0fE2SS/6F8bEc8Aa4HlkiZKWkAyULw+InYAByVdkM4iuoKTxxrMzJpAEJHtaCaNHkB+kYjYJOlW4GHgOHB1RAy8230VcBMwiWTgeFQGj83Mcgvo8wByPmnroPL7amB1les2AOeOUbHMzIYtgChfLGiOYGBm1kqarQsoCweDIZz20lNz5Rs/Z379i6o4TL73Gnr2782VD2Df7nzvGZhZDSUdQHYwMDMrVHhqqZlZu0vGDBwMzMzaW0BfX/lGkB0MzMwK5paBmVmbiwgPIJuZmaeWtpRTXzo9V77xZ56dK9++3nwbaD+281CufAAHnumpf5GZDZtfOjMza3NlXajOwcDMrEgB/V6byMys3QX9HjMwM2tvZX3prCH7GZiZtaxIgkGWowiSPi4pJM2oSFspaYukzZIuynIftwzMzAo2VgPIkuYDbwOerEhbBCwHzgHOBNZJennFvjBVuWVgZlagiKC/rz/TUYAvAJ/g5P3glwFrIqI3IrYCW4Al9W7UFi2DCRJnThxeVaefNSvXs+K0Obny7Tp8LFe+x545mCsfwIEdv8md18yGNoyWwQxJGyq+d0dEd5aMki4Fno6Igb3kB8wD7qv43pOm1dQWwcDMbCxFf+aXSPdExOKhTkpaB1T7C3MV8Eng7dWyVStSvYI4GJiZFSliOMGgzq3iwmrpkl4FLAAGWgVdwEZJS0haApW7bHUB2+s9y8HAzKxAQXHBYMhnRDwI/LYvW9ITwOKI2CNpLfAtSdeRDCAvBNbXu6eDgZlZkQKib3SDQc3HR2ySdCvwMHAcuLreTCJwMDAzK1b003/86Ng+MuLsQd9XA6uHcw8HAzOzgo12N9FoaItgMKFDLJgyflh5pp01O9ez+qfly/fktiO58u3bfThXPoBjh/fnzmtm1Y3FmMFoaItgYGY2ZsItAzMzI+h3MDAza3MFvmcwlhwMzMwKFBH0Hxvb2URFcDAwMyuYWwZmZu3O3URmZoanljav8eM7mDNz8rDynDJ/fv2Lqnhh3JRc+bY9ty9XvpG8Z2BmxUu2vSxkr4Ix1RbBwMxszESM+XIURXAwMDMrUvg9AzOzthc0dtXSvBq2B7KkT0t6WtID6fHOinMrJW2RtFnSRRXp50t6MD33RQ3a683MrOHS2URZjmbS6JbBFyLibyoTJC0ClgPnkGzMsE7Sy9P1uG8AVpDs7/lDYClw59gW2cyslnLOJmpYy6CGZcCaiOiNiK3AFmCJpLnA9Ii4NyIC+AZwWQPLaWZWlVsGw/dRSVcAG4C/iIjngHkkf/kP6EnTjqWfB6dXJWkFSSsCoPftj258aFglu2TjsC4/4aM58w3bDGDPWD1sLDzzwNehBetFa9YJWrNerxjpDeL5PXcd3XjjjIyXN81/v1ENBpLWAXOqnFpF0uXzWZLxls8Cfwt8CKg2DhA10quKiG6gOy3HhohYPKzCN7lWrBO0Zr1asU7QmvWStGGk94iIpUWUZayNajCIiAuzXCfpq8Ad6dceoPKNry5ge5reVSXdzMxGqJGzieZWfH03MNCNsxZYLmmipAXAQmB9ROwADkq6IJ1FdAVw+5gW2sysRTVyzOBzkl5N0tXzBPAnABGxSdKtwMPAceDqdCYRwFXATcAkkllEWWcSdRdW6ubRinWC1qxXK9YJWrNerVinTJRMzDEzs3bWjFNLzcxsjDkYmJlZawcDSUvTJS22SLqm0eUZDklfk7RL0kMVaadLulvSY+nPl1Scq7qERzORNF/Sv0p6RNImSR9L00tbL0mnSFov6ZdpnT6Tppe2TpUkdUr6haQ70u+lr5ekJ9JlbR4YmEraCvUasYhoyQPoBH4D/AdgAvBLYFGjyzWM8r8ReC3wUEXa54Br0s/XAH+dfl6U1m8isCCtd2ej61ClTnOB16afpwG/Tste2nqRvP8yNf08HvgZcEGZ6zSofn8OfAu4oxX+N5iW9QlgxqC00tdrpEcrtwyWAFsi4vGIOAqsIVnqohQi4ifAs4OSlwE3p59v5sRyHFWX8BiLcg5HROyIiI3p54PAIyRvkZe2XpE4lH4dnx5Bies0QFIX8C7gxork0tdrCK1ar8xaORjMA56q+F5z+YqSmB3J+xakP2el6aWrq6SzgdeQ/CVd6nqlXSkPALuAuyOi9HVK/R3wCaBy265WqFcAP5L083TZGmiNeo1Io9cmGk3DWr6i5EpVV0lTge8CfxYRB2qsRF6KekXyHsyrJZ0G3Cbp3BqXl6JOki4GdkXEzyW9OUuWKmlNV6/UGyJiu6RZwN2SHq1xbZnqNSKt3DIYalmLMts58OZ2+nNXml6aukoaTxIIbomI76XJpa8XQETsA35MsrR62ev0BuBSSU+QdLH+Z0n/SPnrRURsT3/uAm4j6fYpfb1GqpWDwf3AQkkLJE0g2SNhbYPLNFJrgSvTz1dyYjmOqkt4NKB8NaXLiPwD8EhEXFdxqrT1kjQzbREgaRJwIfAoJa4TQESsjIiuiDib5P87/xIRH6Dk9ZI0RdK0gc/A20mWwil1vQrR6BHs0TyAd5LMWPkNsKrR5Rlm2b8N7ODE0t0fBs4A7gEeS3+eXnH9qrSem4F3NLr8Q9Tp90ma2L8CHkiPd5a5XsB5wC/SOj0EfCpNL22dqtTxzZyYTVTqepHMLvxlemwa+L1Q9noVcXg5CjMza+luIjMzy8jBwMzMHAzMzMzBwMzMcDAwMzMcDMzMDAcDKxlJn5b08azXSPqgpDOH+YybJF0+knKalY2DgbW6DwLDCgZm7cjBwJqepFXpxiLrgFdUpP+OpP+brj75b5JeOSjf5cBi4JZ0I5NJkj4l6X5JD0nq1tCr5L1R0r9LerxWK0HSVEn3SNqYbphSmmXSzSo5GFhTk3Q+ydo4rwHeA/xuxelu4E8j4nzg48D1lXkj4v8AG4D3R8SrI+IF4MsR8bsRcS4wCbh4iEfPJVk+42Lg2hpFPAK8OyJeC7wF+NsaAcasabXyEtbWGv4AuC0ingeQtDb9ORX4PeCfKn73Tsxwv7dI+gQwGTidZH2aH1S57vsR0Q88LGl2jfsJ+J+S3kiy7v88YDbwTIaymDUNBwMrg2oLaHUA+yLi1VlvIukUktbD4oh4StKngVOGuLy3MmuN274fmAmcHxHH0iWfh7qnWdNyN5E1u58A7077+6cBlwBExAFgq6T3QbI8tqT/VCX/QZL9luHEL+k9acuiiBlDp5JsAnNM0luAlxZwT7Mx55aBNbWI2CjpOyTLXW8D/q3i9PuBGyT9Fcnew2tIliaudBPwFUkvAK8Hvgo8SLIp+v0FFPEW4AeSNqRlrLVrllnT8hLWZmbmbiIzM3M3kVkmkl4FfHNQcm9EvK4R5TErmruJzMzM3URmZuZgYGZmOBiYmRkOBmZmBvx/IFAMz1ezbLwAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# both concentrations on each layer increase\n",
    "phi_a_0 = 500 # age of upper thermocline in CanESM5\n",
    "phi_a_1 = 600\n",
    "\n",
    "phi_b_0 = 2500 # age of lower thermocline in CanESM5\n",
    "phi_b_1 = 2650\n",
    "\n",
    "ha = 500 # thicknes of lower thermocline in CanESM5\n",
    "hb = 2000\n",
    "\n",
    "delta_h_a = np.linspace(0,500, 20)\n",
    "delta_h_b = np.linspace(-500,0, 30)\n",
    "delta_phi_c_results = np.ones([20,30])\n",
    "\n",
    "for ii,dha in enumerate(delta_h_a):\n",
    "    for jj,dhb in enumerate(delta_h_b):\n",
    "        delta_phi_c = phi_c(phi_a_1, phi_b_1, ha+dha, hb+dhb) - phi_c(phi_a_0, phi_b_0, ha, hb)\n",
    "        delta_phi_c_results[ii, jj] = delta_phi_c\n",
    "        \n",
    "plt.pcolormesh(delta_h_a, delta_h_b, delta_phi_c_results.T,cmap='RdBu_r', vmin=-50, vmax=50)\n",
    "plt.colorbar()\n",
    "plt.xlabel('delta h_a')\n",
    "plt.ylabel('delta h_b')\n",
    "plt.axhline(-200)\n",
    "plt.axvline(200)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This actually convinces me. I grabbed the values by eye, but this matches pretty damn good!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python [conda env:conda_tigressdata-omz_paper]",
   "language": "python",
   "name": "conda-env-conda_tigressdata-omz_paper-py"
  },
  "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.8.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}