Parcourir la source

drain tank modeling with pysketcher

Gilbert Brault il y a 5 ans
Parent
commit
d602e264f2
2 fichiers modifiés avec 444 ajouts et 3 suppressions
  1. 441 0
      notebooks/Gravity Drain Tank Modeling.ipynb
  2. 3 3
      notebooks/drainingtank.yml

+ 441 - 0
notebooks/Gravity Drain Tank Modeling.ipynb

@@ -0,0 +1,441 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Gravity Drain Tank Modeling"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 1,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "%matplotlib widget"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 2,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from pysketcher import *"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 3,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from IPython.display import HTML, SVG, display, clear_output"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 4,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from ipywidgets import interact, FloatSlider, Button"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 5,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from math import pi\n",
+    "import numpy as np\n",
+    "from ipywidgets import VBox, FloatSlider, Output, Button, Label"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 6,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import numpy as np\n",
+    "from scipy.integrate import odeint\n",
+    "import time"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Experiment, Simulation (x10)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 7,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "myfig={}\n",
+    "sketch = Sketch(myfig)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 8,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "True"
+      ]
+     },
+     "execution_count": 8,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "sketch.file2Sketch(\"drainingtank.yml\")"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 9,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "def tank(h):\n",
+    "    sketch.container['h'] = h\n",
+    "    sketch.refresh('interieur')\n",
+    "    sketch.refresh('contenu')\n",
+    "    sketch.refresh('hc')\n",
+    "    sketch.refresh('V')\n",
+    "    sketch.refresh('Y')\n",
+    "    sketch.refresh('jet')\n",
+    "    sketch.refresh('dim')\n",
+    "    sketch.refresh('tank')\n",
+    "    drawing_tool.erase()\n",
+    "    sketch.container['tank'].draw()\n",
+    "    display(SVG(Sketch.matplotlib2SVG()))"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 10,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "application/vnd.jupyter.widget-view+json": {
+       "model_id": "5bd9cae3ef5647f98dd8ec3fdb0a5a43",
+       "version_major": 2,
+       "version_minor": 0
+      },
+      "text/plain": [
+       "interactive(children=(FloatSlider(value=1.5, description='h', max=2.0, min=0.12), Output()), _dom_classes=('wi…"
+      ]
+     },
+     "metadata": {},
+     "output_type": "display_data"
+    },
+    {
+     "data": {
+      "text/plain": [
+       "<function __main__.tank(h)>"
+      ]
+     },
+     "metadata": {},
+     "output_type": "display_data"
+    },
+    {
+     "data": {
+      "application/vnd.jupyter.widget-view+json": {
+       "model_id": "960e365964084ff192f76ca86d878111",
+       "version_major": 2,
+       "version_minor": 0
+      },
+      "text/plain": [
+       "Button(description='Simulate', style=ButtonStyle())"
+      ]
+     },
+     "metadata": {},
+     "output_type": "display_data"
+    }
+   ],
+   "source": [
+    "w = FloatSlider(min=sketch.container['d_o']*1.2, max=sketch.container['H'], step=0.1, value=sketch.container['h'])\n",
+    "b = Button(description=\"Simulate\")\n",
+    "i = interact(tank, h=w)\n",
+    "display(i,b)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {
+    "colab_type": "text",
+    "id": "xSIqWFm93OPg",
+    "nbpages": {
+     "level": 2,
+     "link": "[2.2.2 Torricelli's law](https://jckantor.github.io/CBE30338/02.02-Gravity-Drained-Tank.html#2.2.2-Torricelli's-law)",
+     "section": "2.2.2 Torricelli's law"
+    }
+   },
+   "source": [
+    "## Torricelli's law\n",
+    "\n",
+    "Torricelli's law states the velocity of an incompressible liquid stream exiting a liquid tank at level $h$ below the surface is \n",
+    "\n",
+    "$$v = \\sqrt{2gh}$$ \n",
+    "\n",
+    "This is the same velocity as an object dropped from a height $h$. The derivation is straightforward. From Bernoulli's principle,\n",
+    "\n",
+    "$$\\frac{v^2}{2} + gh + \\frac{P}{\\rho} = \\text{constant}$$\n",
+    "\n",
+    "Applying this principle, compare a drop of water just below the surface of the water at distance $h$ above the exit, to a drop of water exiting the tank\n",
+    "\n",
+    "$$gh + \\frac{P_{atm}}{\\rho} = \\frac{v^2}{2} + \\frac{P_{atm}}{\\rho}$$\n",
+    "\n",
+    "$$\\implies v^2 = 2gh$$\n",
+    "$$\\implies v = \\sqrt{2gh}$$\n",
+    "\n",
+    "Torricelli's law is an approximation that doesn't account for the effects of fluid viscosity, the specific flow geometry near the exit, or other flow non-idealities. Nevertheless it is a useful first approximation for flow from a tank."
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {
+    "colab_type": "text",
+    "id": "H8dmoxAJ3OPk",
+    "nbpages": {
+     "level": 2,
+     "link": "[2.2.3 Mass Balance for Tank with Constant Cross-Sectional Area](https://jckantor.github.io/CBE30338/02.02-Gravity-Drained-Tank.html#2.2.3-Mass-Balance-for-Tank-with-Constant-Cross-Sectional-Area)",
+     "section": "2.2.3 Mass Balance for Tank with Constant Cross-Sectional Area"
+    }
+   },
+   "source": [
+    "# Mass Balance for Tank with Constant Cross-Sectional Area\n",
+    "\n",
+    "For a tank with constant cross-sectional area, such as a cylindrical or rectangular tank, the liquid height is described by a differential equation\n",
+    "\n",
+    "$$A\\frac{dh}{dt} = q_{in}(t) - q_{out}(t)$$\n",
+    "\n",
+    "where $$q_{out}$$ is a function of liquid height. Torricelli's law tells the outlet flow from the tank is proportional to square root of the liquid height\n",
+    "\n",
+    "$$ q_{out}(h) = C_v\\sqrt{h} $$\n",
+    "\n",
+    "Dividing by area we obtain a nonlinear ordinary differential equation \n",
+    "\n",
+    "$$ \\frac{dh}{dt} = - \\frac{C_V}{A}\\sqrt{h} + \\frac{1}{A}q_{in}(t) $$\n",
+    "\n",
+    "in our standard form where the LHS derivative appears with a constant coefficient of 1."
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Equation Solving"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Construction parameters"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 11,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "A   = pi * sketch.container['l']**2 / 4              # Tank area [meter^2]\n",
+    "h0 = sketch.container['h']-sketch.container['d_o']   # tank initial height\n",
+    "Cv  = 0.1/60                                         # Outlet valve constant [cubic meters/sec/meter^1/2]"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 12,
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "1.4\n"
+     ]
+    }
+   ],
+   "source": [
+    "print(h0)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Physical Constants"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 13,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "g = sketch.container['g']                         # gravity constant"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Time Horizon"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 14,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "horizon = 10 * 60   # simulation period"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {
+    "colab": {},
+    "colab_type": "code",
+    "id": "qYEqVEsU3OPv",
+    "nbpages": {
+     "level": 3,
+     "link": "[2.2.4.2 Step 2. Define parameters](https://jckantor.github.io/CBE30338/02.02-Gravity-Drained-Tank.html#2.2.4.2-Step-2.-Define-parameters)",
+     "section": "2.2.4.2 Step 2. Define parameters"
+    }
+   },
+   "source": [
+    "## System Solving"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 15,
+   "metadata": {
+    "colab": {},
+    "colab_type": "code",
+    "id": "nuIhX1VO3OPy",
+    "nbpages": {
+     "level": 3,
+     "link": "[2.2.4.3 Step 3. Write Functions for the RHS of the Differential Equations](https://jckantor.github.io/CBE30338/02.02-Gravity-Drained-Tank.html#2.2.4.3-Step-3.-Write-Functions-for-the-RHS-of-the-Differential-Equations)",
+     "section": "2.2.4.3 Step 3. Write Functions for the RHS of the Differential Equations"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "# inlet flow rate in cubic meters/sec\n",
+    "def qin(t):\n",
+    "    return 0.0\n",
+    "\n",
+    "def deriv(h,t):\n",
+    "    return qin(t)/A - Cv*np.sqrt(np.abs(h))/A"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 16,
+   "metadata": {
+    "colab": {},
+    "colab_type": "code",
+    "id": "HiA2HZcb3OP1",
+    "nbpages": {
+     "level": 3,
+     "link": "[2.2.4.4 Step 4. Choose an Initial Condition, a Time Grid, and Integrate the Differential Equation](https://jckantor.github.io/CBE30338/02.02-Gravity-Drained-Tank.html#2.2.4.4-Step-4.-Choose-an-Initial-Condition,-a-Time-Grid,-and-Integrate-the-Differential-Equation)",
+     "section": "2.2.4.4 Step 4. Choose an Initial Condition, a Time Grid, and Integrate the Differential Equation"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "def solver():\n",
+    "    global D,t,h,v0\n",
+    "    D = 2 * np.sqrt(A/pi)\n",
+    "    IC = [h0]\n",
+    "    t = np.linspace(0,horizon,101)\n",
+    "    h = odeint(deriv,IC,t)\n",
+    "    v0 = np.sqrt(2*g*np.abs(h))"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 17,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "solver()"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Animating Gravity Drain Tank"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 18,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "def simulation(b):\n",
+    "    global h, sketch, w\n",
+    "    for hin in [hauteur for hauteur in h if hauteur > sketch.container['d_o']]:\n",
+    "        w.value = hin\n",
+    "        time.sleep(0.1)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 19,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "b.on_click(simulation)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": []
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Python 3",
+   "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.8.2"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 4
+}

+ 3 - 3
notebooks/drainingtank.yml

@@ -16,15 +16,15 @@
         g: 9.81       # gravity constant 
   - name: variables
     shapes:
-        h: 1.0        # Current tank heigth
+        h: 1.5        # Current tank heigth
         V: sqrt(2*g*(h-d_o))
-        X: np.linspace(0.0, 2.0, 20)
+        X: np.linspace(0.0, 1.2, 12)
         Y: -g/2*(X/V)**2
   - name: frame
     shapes:
         setframe:
             action: drawing_tool.set_coordinate_system(xmin=-l*1.2, xmax=l*3,
-                                   ymin=-H, ymax= H,
+                                   ymin=-2.8, ymax= H,
                                    axis=False)
   - name: scene
     shapes: