{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "3a8da1a3-5b0d-440c-b1d5-ab676a8a0e4b",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pylab as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "871068ea-2174-432c-b46b-257d6ae114ef",
   "metadata": {},
   "source": [
    "### Solving the expression $1 - 2e^{-(R/h)} + R/h = 0.5$\n",
    "First, let's simplify things with the substitution $x=R/h$ so that the equation \n",
    "becomes $1-2e^{-x}+x=0.5$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "c9846a25-55aa-42f7-97c7-c16dab2a5094",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0, 0.5, 'f(x)')"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAEGCAYAAABsLkJ6AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAhzElEQVR4nO3dd3xV5eHH8c9DCIS9ZyAJEHaQFYZQF+5V1GqjKKgoqHVba22rta3Van9WK7ZqWGUooNZRUaRWlGkYCXsEQhICIUAWZK+bPL8/kqaUglzJOPfmfN+vF68XJ/cm93s44X7vOec5zzHWWkRExH0aOR1AREScoQIQEXEpFYCIiEupAEREXEoFICLiUo2dDvB9dOzY0YaFhTkdQ0TEr8TFxWVaazud+nW/KoCwsDBiY2OdjiEi4leMMSmn+7oOAYmIuJQKQETEpVQAIiIupQIQEXEpFYCIiEupAEREXEoFICLiUioAEREflpJVwO+W7uZEYWmt/2y/uhBMRMQNrLVsSM5mztpkvtpzjABjGNenA5cN6lKrr6MCEBHxEaWeCj7bnsactcnsSsulXfNAHrw4nMnnh9KldVCtv54KQETEYdkFpby7PoUF61PIyCshvHNLXrxxCDcOD6ZZk4A6e10VgIiIQxKO5TF3XTIfbT5MiaeCC/t14pVbenFh344YY+r89VUAIiL1yFrL6oRM5qxNZvW+DJo2bsRNI3owdXwYfbu0qtcsKgARkXpQXFbOR5sPM3ddMvvT8+ncqilPXtGPSWNCad+iiSOZVAAiInUoPbeYBTEpvLshheOFZUQEt+a1qKFcO6Q7TRo7OxJfBSAiUgfij+Yya3Uyn247jKfCcvnALtzzg16M7tW+Xo7ve0MFICJSS/59fH/2miTWJGTSLDCASaNDmPqDXoR2aOF0vP+hAhARqaESTzn/2JrGnDXJ7D2WR+dWTfnZlf25fUwIbZs7c3zfGyoAEZFzdKKwlHc3HGTetwfIyCthQNdWvHLLUH441Pnj+95QAYiIfE8HMguYuy6ZD2JTKSor58J+nXj1x734QXj9jN+vLSoAEREvWGuJSznOrDVJfLn7GIGNGjFxWHfuvaA3/bvW7/j92uJYARhjgoDVQNOqHH+31j7nVB4RkdPxlFfwz13HmLUmia2HTtC2an6eKeNC6dyq9ufnqU9O7gGUABOstfnGmEBgrTHmC2vtegcziYgAkF/i4f1Nh5i7LpnU40WEdWjO8xMH86ORPWjepGEcPHFsLay1FsivWgys+mOdyiMiAnAkp4h56w6waONB8oo9RIa249nrBnHZwC4ENPKf4/vecLTGjDEBQBwQDvzVWrvhNM+ZDkwHCAkJqd+AIuIae47kMnN1Eku3pVFhLVdHdOPeC3oxPKSd09HqjKMFYK0tB4YZY9oCHxtjIqy1O095zkxgJkBkZKT2EESk1lhriUnKInpVEqv2ZdC8SQCTzw9l6vhe9Gzf3Ol4dc4nDmRZa08YY1YCVwE7z/J0EZEaKa+wLN95lOjViWxPzaFjyyY8eUU/7hgb6tMXbtU2J0cBdQLKqt78mwGXAS87lUdEGr7isnI+iEtl9pokUrIKCevQnBdvHMJNI4IJCqy7G6/4Kif3ALoB86vOAzQC3rfWfuZgHhFpoI4XlLJwfQrzvz1AVkEpQ3u25RdXD+DyQV0b3Ind78PJUUDbgeFOvb6INHyHsguZszaZ9zYdoqisnAkDOnPfhb19akZOJ/nEOQARkdq0Ky2HmauT+Gz7EQwwcVgw0y/03yt264oKQEQaBGst3yZm8faqRNYkZNKiSQBTx4dx9/hedG/bzOl4PkkFICJ+zVNewRdVI3p2Hs6lY8umPHVVf24fE0qbZoFOx/NpKgAR8Uv/HtEzc3Uih7KL6N2xBS/dNIQbhrtzRM+5UAGIiF/JKy7j3Q0Hmb0mmcz8Eob1bMsz1w7i8oFdaOTiET3nQgUgIn4hK7+Ev607wIKYA+QWe7igb0ceuHgY5/fuoBE950gFICI+7fCJImatTmLJpoOUeCq4anBXHri4D+f1aOt0NL+nAhARn7Q/PZ+3VyXyyZbDANwwPJj7L+pDeOeWDidrOFQAIuJTdqTm8ObK/SzfdZSmjRtxx9hQpl3Ym2AN5ax1KgARcZy1lvVJ2by5cj9rEjJpFdSYBy8O5+7xYXRo2dTpeA2WCkBEHFNRYVkRn86bK/ez5eAJOrZsys+vGsAdY0NoFaQx/HVNBSAi9c5TXsFn24/w5sr97DuWT492zXj+hghuGdlDY/jrkQpAROpNqaeCjzan8ubKRA5mF9KvS0teixrK9ed1p3FAI6fjuY4KQETqXHFZOR/EHuLtVUkcPlHEeT3a8KtrR+riLYepAESkzhSVlrNo40GiVyWSnlfCyNB2vHBjBBf166SLt3yACkBEal1+iYd31qcwe00SmfmljO3dnj9HDeP8Prpq15eoAESk1uQUlTH/2wPMXZfMicIyLujbkUcu7cuosPZOR5PTUAGISI0dLyhl7rpk5q07QF6Jh8sGduahCX0Z1rOt09HkO6gAROScZeaXMGtNEu/EpFBQWs7VEV15aEI4g7u3cTqaeEEFICLf27HcYqJXJbFoYwqlngquO687D00Ip18X3XLRn6gARMRr6bnFvLkykUUbD1JeYblhWDAPXtKH3p00QZs/UgGIyFml5xbz1qpEFm04iKfCcvOIHjx4STghHZo7HU1qQAUgImeUkVfC26sSeWd9Cp4Ky03Dg3l4Ql+98TcQKgAR+R+Z+SVEr0pk4frKY/w3Du/BwxPCCevYwuloUotUACJSLTO/hJmrk1gYk0KJp5wbqj7x99Ibf4OkAhARsvJLmLkmiQXfVr7xTxwWzMMTwnVyt4FTAYi4WHZBKTNXJ7Eg5gBFZeVMHNqdhy/tSx+98buCCkDEhXIKy5i5JpF56w5QWFbO9ed155FL++p+uy6jAhBxkfwSD3PXJjNrTRL5JR6uHdKNRy/tS19dwOVKKgARFyguK2dhTApvrUoku6CUywd14adX9GNA19ZORxMHqQBEGrBSTwXvxx7ija8TOJZbwgV9O/LTK/prkjYBHCwAY0xPYAHQFagAZlprX3cqj0hDUl5h+WTLYf68Yh+HsouIDG3H67cOZ2zvDk5HEx/i5B6AB/iptXazMaYVEGeM+Ze1dreDmUT8WkWF5YudR3n1X3tJzCggIrg1v7s7got1By45DccKwFp7BDhS9fc8Y8weIBhQAYh8T9Zavtmbzp++3MeutFzCO7fkrdtHcFVEV73xyxn5xDkAY0wYMBzYcJrHpgPTAUJCQuo3mIgfiEnM4pUv9xKXcpyQ9s159cdDmTgsmADdbF3OwvECMMa0BD4EHrPW5p76uLV2JjATIDIy0tZzPBGftfNwDi8vj2dNQiZdWwfxwo0R/DiyJ4EBjZyOJn7C0QIwxgRS+eb/rrX2IyeziPiLg1mFvPLlXj7dlkbb5oH86pqBTD4/lKDAAKejiZ9xchSQAeYAe6y1rzqVQ8RfZOaX8Jev9/PuhhQCGhkevKQP0y/sQ5tmgU5HEz/l5B7AeGAysMMYs7Xqa7+01i5zLpKI78kv8TB7TRKzVidR7KkgalRPHr20L11aBzkdTfyck6OA1gI6SyVyBqWeChZvPMiMFQlkFZRydURXnryyvyZqk1rj+ElgEflvFRWWpdvT+NOX+ziYXcjY3u2ZfdUAhoe0czqaNDAqABEfsiYhg5e+iGdXWi4DurZi3t2juEgXcUkdUQGI+IDdabn84Ys9rEnIpEe7ZrwWNZSJQ4NppLH8UodUACIOOppTzJ++3MvfN6fSplkgz143iDvGhtC0sYZ0St1TAYg4oKDEQ/TqypE95RWWaRf05sGLw2nTXEM6pf6oAETqUXmF5f3YQ7z6r31k5JVw3XndeOrKAYR0aO50NHEhFYBIPVm5N50/LItn77E8Roa2I3rySEZoZI84SAUgUsf2HMnlxWWVJ3hDOzTXLJ3iM1QAInXkWG7lCd4P4lJpHVR5gnfy2FCaNNZkbeIbVAAitay4rJzZa5L46zeJeCoquGd8Lx6aEE7b5k2cjibyX1QAIrXE2sq7cb24bA+px4u4cnAXfnnNQEI7tHA6mshpqQBEasGutBx+t3Q3G5KzGdC1FYvuHcO48I5OxxL5TioAkRrIyi/hlS/3sWTTQdo2C+T5GyK4bVRPGuumLOIHVAAi56DUU8GCmAO8viKBotJy7h7Xi0cv7asLucSvqABEvqdv4tN5/rPdJGUWcFG/Tjx73SDCO2uKZvE/KgARL+1Pz+f5z3azal8GvTu14G93jeKSAZ2djiVyzlQAImeRX+JhxooE5q5NplmTAJ65diBTzg/TeH7xeyoAkTOw1vLptjReXLaHY7klREX25Kmr+tOhZVOno4nUChWAyGnsPZrHr/+xkw3J2QwJbsPbd4zUHbmkwVEBiJwkr7iMP3+VwLxvD9AqqDEv3BjBraNCCNCNWaQBUgGIUHm455Oth3lxWTyZ+SXcOiqEp67sT7sWmr5BGi4VgLjeniO5PPePXWw8kM3QHm2YPSWSoT3bOh1LpM6pAMS1covLeO1f+1gQk0LroMa8dNMQfhzZU/fhFddQAYjrWGtZtuMov126i4z8EiaNDuFnV/bXbJ3iOioAcZVD2YU8+4+drNybweDurZmlwz3iYioAcYWy8gpmr0nm9RX7aGQMz1w7kLvGhWnSNnE1FYA0eHEp2fzyo53sPZbHFYO68JsfDqZ722ZOxxJxnApAGqwThaW8vDyexRsP0b1NELOmRHL5oC5OxxLxGSoAaXD+Pab/95/t4URRGdMu6MVjl/WjRVP9uoucTP8jpEE5kFnArz7Zwbr9WQzt2ZYFN0YwuHsbp2OJ+CQVgDQInvIK5qxN5tV/7aNJQCOenziYSWNCNYWDyHdQAYjf252Wy88/3M6OwzlcPqgLz0+MoGubIKdjifg8rwrAGNMZGA90B4qAnUCstbaiJi9ujJkLXAekW2sjavKzxH2Ky8p54+sEolcl0bZ5IH+dNIJrhnTFGH3qF/HGdxaAMeYS4GmgPbAFSAeCgBuAPsaYvwN/stbmnuPrzwP+Aiw4x+8Xl9p0IJuff7idpIwCfjSiB89cO1ATt4l8T2fbA7gGmGatPXjqA8aYxlR+er8c+PBcXtxau9oYE3Yu3yvulFdcxh+X72Xh+hSC2zZj/tTRXNSvk9OxRPzSdxaAtfZn3/GYB/iktgOdyhgzHZgOEBISUtcvJz7sm/h0fvXxDo7kFnP3+DCevKK/hnaK1IBX18EbYxYaY9qctBxmjFlRd7H+w1o701obaa2N7NRJn/TcKKewjCfe28rd8zbRomlj/n7/OJ67frDe/EVqyNv/QWuBDcaYJ4Bg4GfAT+sslUiVr+OP8fSHO8gqKOWRCeE8OCGcpo0DnI4l0iB4VQDW2mhjzC7gGyATGG6tPVqnycTVcovLeH7pbj6IS6V/l1bMvWsUEcG6oEukNnk7DHQy8CwwBTgPWGaMudtau60mL26MWQxcDHQ0xqQCz1lr59TkZ4r/W70vg59/uJ1jucU8eEkfHrm0rz71i9QBbw8B/Qj4gbU2HVhsjPmYyiGcw2vy4tba22ry/dKw5Jd4eOHzPSzeeJDwzi356CfjGaa5+kXqjLeHgG44ZXmjMWZMnSQSV1q3P5On/r6dIzlF3HdRbx6/rB9BgfrUL1KXznYh2DPAm9ba7FMfs9aWGmMmAM2ttZ/VVUBp2ApKPLz0RTwL16fQq2MLPrh/HCND2zkdS8QVzrYHsANYaowpBjYDGVReCdwXGAZ8BbxYlwGl4YpLOc4T72/lYHYh9/ygF09e0Z9mTfSpX6S+nK0AbrbWjjfGPEXlNBDdgFzgHWC6tbaorgNKw1NWXsEbKxL4yzf76damGUumjWVM7w5OxxJxnbMVwEhjTChwO3DJKY81o3JiOBGvJWbk88R7W9mWmsOPRvTgNz8cRKugQKdjibjS2QrgbWA50BuIPenrBrBVXxc5K2st72w4yAuf7yYoMIC3bh/B1UO6OR1LxNXONhfQDGCGMeYta+0D9ZRJGpj0vGKe/nAHX8enc0Hfjrxyy1C6tNZ8/SJO83YYqN785Zx8uesoT3+0g4ISD7+5fhBTzg+jke7SJeITNJuW1InCUg+/W7qbJZsOMbh7a/4cNYy+XVo5HUtETqICkFq3Ky2HhxdvITmzgAcu7sPjl/WjSWOvJp4VkXqkApBaY61l3rcH+MOyeNo2D+Tde8YwLryj07FE5AxUAFIrsgtKeerv2/hqTzqXDujM/90ylPa6RaOIT1MBSI19m5jJ4+9t5XhBGc9dP4i7xoXpxuwifkAFIOfMU17B61VX9Pbq0II5d2rOfhF/ogKQc5J6vJBHl2wlLuU4t4zswW9+qFs0ivgb/Y+V7+3LXUd58oNtVFh4/dZhTBwW7HQkETkHKgDxWll5BX9cHs+sNckMCW7DXyeNIKRDc6djicg5UgGIV47kFPHQoi3EpRxn8thQnrluoG7TKOLnVAByVqv3ZfDYe1spKStnxm3D+eHQ7k5HEpFaoAKQMyqvsLz+1T7e+GY//Tq34s07RtCnU0unY4lILVEByGll5JXw6JItfJuYxc0je/D8xAjdrUukgVEByP/YkJTFw4u3kFNUxh9vPo8fR/Z0OpKI1AEVgFSz1jJ7TTIvLY8ntH1z5k8dzcBurZ2OJSJ1RAUgABSUeHjqw+18vv0IV0d05Y83n6dbNYo0cCoA4UBmAfctjCMhPY+fXzWA+y/qrbl8RFxABeByX8cf49ElWwloZJg/dTQX9O3kdCQRqScqAJeqqLC88fV+/rxiH4O6tebtO0bSs72u6hVxExWAC+UWl/HEe1v5ak86Nw0P5sWbhhAUqCGeIm6jAnCZhGN5TF8Yx6HsQn77w8FMOT9Ux/tFXEoF4CLfxKfz8OItBAUGsGjaWEb3au90JBFxkArABay1zFydxEvL4xncvTUzJ0fSvW0zp2OJiMMaOfnixpirjDF7jTH7jTFPO5mloSouK+en72/jD1/Ec82Qbnxw3zi9+YsI4OAegDEmAPgrcDmQCmwyxnxqrd3tVKaGJj23mPveiWPLwRM8cXk/Hp4QruP9IlLNyT2A0cB+a22StbYUWAJM/K5vSMoo4IPYQ0DlzUmiomP4eEsqAEWl5URFx7B0WxpQOdIlKjqG5TuPAJBdUEpUdAxf7T4GQHpeMVHRMazcmw5A2okioqJjWJuQCcDBrEKiomNYn5QFQGJGPlHRMcSlZAOw92geUdExbDt0AoBdaTlERcewKy0HgG2HThAVHcPeo3kAxKVkExUdQ2JGPgDrk7KIio7hYFYhAGsTMomKjiHtRBEAK/emExUdQ3peMQBf7T5GVHQM2QWlACzfeYSo6Bhyi8sAWLotjajoGIpKywGYsSKB8S9/zZ60XN66fQTd2gRx68z11f+Wizce5PbZ/1leGHOAO+durF6euzaZe+dvql6euTqR+xfGVS+/uXI/Dy3aXL08Y0UCjy3ZUr386pd7efKDbdXLLy+P5xcfba9efuHz3Tz7yc7q5d8u3cVvl+6qXn72k5288Pl/Pgv84qPtvLw8vnr5yQ+28eqXe6uXH1uyhRkrEqqXH1q0mTdX7q9evn9hHDNXJ1Yv3zt/E3PXJlcv3zl3IwtjDlQv3z57PYs3HqxejoqO0e+el797H29JJSo6hrLyCgA+iD1EVHRM9b+lfvfq/3fvTJwsgGDg0EnLqVVf+y/GmOnGmFhjTGxZWVm9hfNnS7el8cbXCRgM7903lquHdHM6koj4IGOtdeaFjbkFuNJae2/V8mRgtLX24TN9T2RkpI2Nja2viH7HWsufv0rg9RUJRIa24+3JI+nYsqnTsUTEYcaYOGtt5Klfd3IUUCpw8jzDPYA0h7L4vRJPOU9/uIOPtxzm5pE9eOHGCN2yUUS+k5MFsAnoa4zpBRwGbgUmOZjHbx0vKOW+hXFsPJDNz67sz08u7qOTvSJyVo4VgLXWY4x5CPgnEADMtdbuOsu3ySkOZBZw97xNHD5RpPv1isj34uiFYNbaZcAyJzP4s9gD2UxbUHlOZNG9Y4gM05W9IuI9XQnspz7dlsaT728juF0z/nbXKMI6tnA6koj4GRWAn7HW8tdv9vPKl/sY3as90XeMpF2LJk7HEhE/pALwI57yCn758Q7ej03lxuHBvPSjIRrpIyLnTAXgJwpLPTy0aAtfx6fzyIRwHr+8n0b6iEiNqAD8QHZBKVPnbWJ76gleuDGC28eEOh1JRBoAFYCPO5RdyJ1zN3L4RBFv3TGSKwd3dTqSiDQQKgAftjstlzv/tpGSsnLeuXcMozTMU0RqkQrAR32bmMl9C+JoGdSYdx8YR78urZyOJCINjArAB322PY0n3ttGaIfmzJ86WjdwEZE6oQLwMfPWJfPbz3YTGdqO2VNG0aZ5oNORRKSBUgH4CGstM1bs57Wv9nHFoC7MuG04QYEa4y8idUcF4AOstfz+8z3MWZvMzSN78NJNQ2gc4OjtmkXEBVQADiuvsPzio+28H5vKXePC+PV1g2jUSBd4iUjdUwE4qMRTzuPvbWXZjqM8cmlfHr+sr67uFZF6owJwSGGph/vf2czqfRk8c+1A7r2gt9ORRMRlVAAOyCkq4555m9h88Dgv/2gIUaNCnI4kIi6kAqhn2QWlTJ6zgX3H8njjthFce143pyOJiEupAOpRZn4Jd8zeQHJmATOnRHJJ/85ORxIRF1MB1JP03GImzd5A6vFC5t41ivHhHZ2OJCIupwKoB0dzipk0az1Hc4uZd/doxvbu4HQkEREVQF07fKKISbPWk5VfyoKpo3XjdhHxGSqAOnQou5BbZ64nt7iMd+4dw7CebZ2OJCJSTQVQRw5kFnDbrPUUlZWzeNpYIoLbOB1JROS/qADqQHJmAbfOjKGs3LLo3rEM6t7a6UgiIv9DBVDLUrIKuG3mejzllsXTxtK/q27kIiK+SVNO1qJD2YVMmrWBEk85704bozd/EfFp2gOoJYdPFHHbrPXkl3hYNG0MA7rqsI+I+DbtAdSCIzmVQz1zisp4554xDO6uE74i4vtUADV0LLeYSbM2VI/zH9JDb/4i4h9UADWQkVfCpFnrSc8tZv7UUQwPaed0JBERr+kcwDnKLijl9tnrSTtRzPypoxkZqit8RcS/aA/gHOQWlzFl7gZSsgqZc1cko3vpzV9E/I8jBWCMucUYs8sYU2GMiXQiw7kqKi3n3nmxxB/J4+07RjKuj2b1FBH/5NQewE7gJmC1Q69/Tko9Fdz/ThybUrJ5LWoYlwzQfP4i4r8cOQdgrd0D+NUN0MsrLI+/t5VV+zJ46aYhXD+0u9ORRERqxOfPARhjphtjYo0xsRkZGY5ksNbyi4+28/mOI/zqmoHcOlr38BUR/1dnewDGmK+Arqd56FfW2n94+3OstTOBmQCRkZG2luJ5zVrL7z/fw/uxqTwyIZxpF/au7wgiInWizgrAWntZXf3s+jRjxX7mrE3mrnFhPH55P6fjiIjUGp8/BOSkBTEHeO2rfdw8sge/vm6QX52zEBE5G6eGgd5ojEkFzgc+N8b804kc32XZjiM89+kuLhvYhZduGkKjRnrzF5GGxalRQB8DHzvx2t6ISczisSVbGRHSjr9MGk7jAO0oiUjDo3e2U+w5ksv0BbGEdGjOnDsjCQoMcDqSiEidUAGcJPV4IXfO3UiLpo1ZMHU0bZs3cTqSiEidUQFUOV5QypS5GykuK2fBPaPp3raZ05FEROqUZgOlcn6fqfM3kXq8iHfuGUO/LrqVo4g0fK7fA/CUV/DQos1sO3SCGbcO18yeIuIart4DsNby6093sSI+nd/fEMFVEae7cFlEpGFy9R7ArDVJLNpwkAcu7sMdY0OdjiMiUq9cWwBf7DjCi8viufa8bvzsiv5OxxERqXeuLIAtB4/z2HtbGRHSlj/dMlRX+YqIK7muAA5lFzJtQSxdWgcxa4ou9BIR93LVSeCcojKmzttEqaeCJdNH0aFlU6cjiYg4xjUFUFZewU/ejSM5s4AF94wmvHNLpyOJiDjKFQVgreWZj3eybn8W/3fzebqRu4gILjkH8PaqJN6LPcTDE8K5JbKn03FERHyCKwqgR7tm3DKyB0/ojl4iItVccQjo+qHduX5od6djiIj4FFfsAYiIyP9SAYiIuJQKQETEpVQAIiIupQIQEXEpFYCIiEupAEREXEoFICLiUsZa63QGrxljMoCUc/z2jkBmLcZxktbFN2ldfJPWBUKttZ1O/aJfFUBNGGNirbWRTueoDVoX36R18U1alzPTISAREZdSAYiIuJSbCmCm0wFqkdbFN2ldfJPW5Qxccw5ARET+m5v2AERE5CQqABERl2pwBWCMucoYs9cYs98Y8/RpHjfGmBlVj283xoxwIqc3vFiXi40xOcaYrVV/fu1EzrMxxsw1xqQbY3ae4XF/2iZnWxe/2CYAxpiexphvjDF7jDG7jDGPnuY5Pr9tvFwPv9guxpggY8xGY8y2qnX57WmeU3vbxFrbYP4AAUAi0BtoAmwDBp3ynGuALwADjAU2OJ27ButyMfCZ01m9WJcLgRHAzjM87hfbxMt18YttUpW1GzCi6u+tgH3++P/Fy/Xwi+1S9e/csurvgcAGYGxdbZOGtgcwGthvrU2y1pYCS4CJpzxnIrDAVloPtDXGdKvvoF7wZl38grV2NZD9HU/xl23izbr4DWvtEWvt5qq/5wF7gOBTnubz28bL9fALVf/O+VWLgVV/Th2pU2vbpKEVQDBw6KTlVP73F8Gb5/gCb3OeX7W7+IUxZnD9RKt1/rJNvOV328QYEwYMp/IT58n8att8x3qAn2wXY0yAMWYrkA78y1pbZ9ukod0U3pzma6e2pzfP8QXe5NxM5Rwf+caYa4BPgL51HawO+Ms28YbfbRNjTEvgQ+Axa23uqQ+f5lt8ctucZT38ZrtYa8uBYcaYtsDHxpgIa+3J55xqbZs0tD2AVKDnScs9gLRzeI4vOGtOa23uv3cXrbXLgEBjTMf6i1hr/GWbnJW/bRNjTCCVb5rvWms/Os1T/GLbnG09/G27AFhrTwArgatOeajWtklDK4BNQF9jTC9jTBPgVuDTU57zKTCl6kz6WCDHWnukvoN64azrYozpaowxVX8fTeX2zKr3pDXnL9vkrPxpm1TlnAPssda+eoan+fy28WY9/GW7GGM6VX3yxxjTDLgMiD/labW2TRrUISBrrccY8xDwTypH0cy11u4yxtxf9fjbwDIqz6LvBwqBu53K+128XJebgQeMMR6gCLjVVg0T8CXGmMVUjsLoaIxJBZ6j8uSWX20T8Gpd/GKbVBkPTAZ2VB1zBvglEAJ+tW28WQ9/2S7dgPnGmAAqS+p9a+1ndfUepqkgRERcqqEdAhIRES+pAEREXEoFICLiUioAERGXUgGIiLiUCkBExKVUACIiLqUCEKkBY8yoqjnZg4wxLarmcI9wOpeIN3QhmEgNGWN+DwQBzYBUa+0fHI4k4hUVgEgNVc3VtAkoBsZVzeYo4vN0CEik5toDLam8G1WQw1lEvKY9AJEaMsZ8SuUd23oB3ay1DzkcScQrDWo2UJH6ZoyZAnistYuqZnD81hgzwVr7tdPZRM5GewAiIi6lcwAiIi6lAhARcSkVgIiIS6kARERcSgUgIuJSKgAREZdSAYiIuNT/A+T6ueciDAieAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# note that 1 - 2e^(-x) + x = 0.5 is the same as 1 - 2e^(-x) + x - 0.5 = 0\n",
    "# so we want to know the value of x where  1 - 2e^(-x) + x - 0.5 = 0.\n",
    "# this is called \"root-finding\".\n",
    "\n",
    "# define that function\n",
    "def func(x):\n",
    "    return 1 - 2*np.exp(-x) + x - 0.5\n",
    "\n",
    "# and make a quick plot\n",
    "x = np.arange(0,3,0.01)\n",
    "plt.plot(x,func(x))\n",
    "plt.gca().axhline(y=0,ls=':')\n",
    "plt.xlabel('x')\n",
    "plt.ylabel('f(x)')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d4c3a6f8-bfb7-40e4-b2de-d4736338fbf0",
   "metadata": {},
   "source": [
    "### root-finding by determining the value in the array f(x) which is closest to 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "b2918783-226e-4bc9-bb1c-5b0de84e1ef9",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Solution is x = 0.599\n"
     ]
    }
   ],
   "source": [
    "# evaluate t over a range of values that you think bracket the solution\n",
    "# the step-size will determine the accuracy of your answer\n",
    "x = np.arange(0,3,0.001) \n",
    "\n",
    "# find the index in the t array where func is closest to 0.\n",
    "idx = np.argmin(np.abs(func(x))) \n",
    "\n",
    "# x[idx] is your solution\n",
    "solution1=x[idx]\n",
    "print('Solution is x = {:.3f}'.format(solution1)) "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "371c2389-0f9b-4383-a1e3-9219239f4b82",
   "metadata": {},
   "source": [
    "### root-finding by using scipy's root-finding routines"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "c063ae10-a840-4a01-abc6-9a59ded00eda",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Printing out all the info from the root-finder:\n",
      "*****************************************\n",
      "    fjac: array([[-1.]])\n",
      "     fun: array([-2.22044605e-16])\n",
      " message: 'The solution converged.'\n",
      "    nfev: 7\n",
      "     qtf: array([1.70319314e-12])\n",
      "       r: array([-2.09886731])\n",
      "  status: 1\n",
      " success: True\n",
      "       x: array([0.59886728])\n",
      "*****************************************\n",
      "Solution is x = 0.599\n"
     ]
    }
   ],
   "source": [
    "from scipy.optimize import root\n",
    "\n",
    "# you need to give the root-finder a rough guess at the correct value\n",
    "sol_guess = 0.5  \n",
    "\n",
    "# call the root-finding routine, printing out all its info\n",
    "print('Printing out all the info from the root-finder:')\n",
    "print('*****************************************')\n",
    "print(root(func,sol_guess)) # give it the function and your initial guess for the solution\n",
    "print('*****************************************')\n",
    "\n",
    "# the solution is given in the x parameter of the root-finder, and we will\n",
    "# take the first root value (which in this case is the only root value\n",
    "solution2 = root(func,sol_guess).x[0]\n",
    "print('Solution is x = {:.3f}'.format(solution2)) "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a1c7f257-0fe9-4398-8fca-0e809022c8e3",
   "metadata": {},
   "source": [
    "### So both ways give a solution of x=0.599. \n",
    "### Since $x=R/h$, this means that the solution is $R = 0.599h$"
   ]
  }
 ],
 "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.9.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}