{ "metadata": { "name": "", "signature": "sha256:b7fad7b1addf530baa9cb1abca4446bf0d616d304170a5a6cb62ed37999f3ba8" }, "nbformat": 3, "nbformat_minor": 0, "worksheets": [ { "cells": [ { "cell_type": "code", "collapsed": false, "input": [ "import numpy as np\n", "import matplotlib.pyplot as plt\n", "import bisect\n", "from IPython.html.widgets import interact\n", "\n", "# example based on http://lorenabarba.com/blog/cfd-python-12-steps-to-navier-stokes/" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 1 }, { "cell_type": "code", "collapsed": false, "input": [ "c = 1.0" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 2 }, { "cell_type": "code", "collapsed": false, "input": [ "# example showing numerical diffusion and osciliations\n", "def plot(nx, dt, nt, bc):\n", " # create the x grid\n", " x = np.arange(0,1,1.0/nx)\n", " # equidistant grid\n", " dx = x[1] - x[0]\n", " # initial waterlevel\n", " s1 = np.zeros(nx) + 0.1\n", " # use current waterlevel as previous waterlevel\n", " s0 = s1.copy()\n", " for n in range(nt): \n", " # set the boundary condition\n", " s1[0] = bc\n", " # upwind scheme\n", " s1[1:] = s0[1:] - c*dt/dx*(s0[1:] - s0[:-1])\n", " # set previous to current\n", " s0 = s1.copy()\n", " plt.step(x, s1)\n", " plt.ylim(0, 3)\n", "_ = interact(plot, nx=(10, 100), dt=(0.005, 0.02,0.001), bc=(0, 2.0), nt=(0,50))" ], "language": "python", "metadata": {}, "outputs": [ { "metadata": {}, "output_type": "display_data", "png": "iVBORw0KGgoAAAANSUhEUgAAAXgAAAECCAYAAAD0JMwBAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAFZdJREFUeJzt3X9MVff9x/HXBdKLIvbu3qW6SRysZamamal0NIHNtQbj\nuiypY/7RLqZuNV1wU7AbZdptoX84A05rXQHbKGFp9kf/aLr0jy2W7EfqbWcWxiV1zK9XVtvMtIPA\n7S2VqS33nu8fnfcjAy73crj33Ht4Pv7ywzlc3r7Edw7v+zkHj2VZlgAArlPgdAEAgMygwQOAS9Hg\nAcClaPAA4FI0eABwKRo8ALhUUbKDH330kVpbW/Xxxx9rcnJS99xzjx5++OFp53V3d2tgYEBer1d7\n9uxRRUVFxgoGAKTGM9c++Bs3bsjr9SoWi+nnP/+5du7cqbvvvjtxvL+/X2fOnNGBAwd06dIl9fT0\n6NChQxkvHACQ3JwjGq/XK0manJxUPB7XsmXLphzv6+vT5s2bJUmVlZWamJhQNBrNQKkAgHQkHdFI\nUjweV0tLi4aHh7V161aVlZVNOR6JRBQIBBLrQCCgSCQin8+38NUCAFI25xV8QUGBjhw5opMnT+rC\nhQsaHBycdg5POwCA3DPnFfxNS5cu1YYNG/TPf/5T69atS3zc7/drbGwssR4bG5Pf75/1dV599VUV\nFhbOs1wAWJx8Pp82bdqU1uckbfDj4+MqLCxUSUmJPvroI50/f17f/va3p5xTVVWlM2fOqKamRuFw\nWCUlJUnHM4WFhdq4cWNaRQLAYtff35/25yRt8NFoVB0dHYrH47IsS1/96lf1xS9+Ub29vZKkuro6\nbdy4UaFQSHv37lVxcbEaGhrmV/0iFAwGVVtb63QZOYEsDLIwyMKepA1+9erVamtrm/bxurq6KetH\nH310YasCANg25z74hfaHP/yBEQ0ApKm/v19btmxJ63N4VAEAuBQN3kHBYNDpEnIGWRhkYZCFPTR4\nAHApZvAAkAeYwQMAEmjwDmK+aJCFQRYGWdhDgwcAl2IGDwB5gBk8ACCBBu8g5osGWRhkYZCFPTR4\nAHApZvAAkAeYwQMAEmjwDmK+aJCFQRYGWdhDgwcAl2IGDwB5gBk8ACCBBu8g5osGWRhkYZCFPTR4\nAHApZvAAkAeYwQMAEmjwDmK+aJCFQRYGWdhDgwcAl2IGDwB5gBk8ACCBBu8g5osGWRhkYZCFPTR4\nAHCppDP40dFRdXR06IMPPpDH49GWLVv0wAMPTDlncHBQ7e3tWrFihSSpurpa9fX1s35BZvAAkL75\nzOCLkh4sKtIjjzyi8vJyXb9+XS0tLVq/fr3KysqmnLd27Vq1tLSkXzEAIGOSjmh8Pp/Ky8slScXF\nxVq1apXef//9aedleSOOazBfNMjCIAuDLOxJegV/q5GREb399tuqrKyc8nGPx6NwOKzm5mb5/X7t\n3Llz2hU+ACD7UtoHf/36dbW2tupb3/qWvvzlL085du3aNRUUFMjr9SoUCqmnp0fPPPPMrK/FDB4A\n0peRffCTk5M6evSovvKVr0xr7pK0ZMkSeb1eSdKGDRs0OTmpq1evJn3NW3/sCgaDrFmzZs06hXW6\nkl7BW5aljo4OLVu2TLt27ZrxnGg0qttvv10ej0dDQ0N6+umn1dHRMesX5AreCAaDqq2tdbqMnEAW\nBlkYZGEs+C6aixcv6uzZs1q9erWeeOIJSdJDDz2k0dFRSVJdXZ3OnTun3t7exJimsbFxnuUDABYS\nz6IBgDzAs2gAAAk0eAfZefPEbcjCIAuDLOyhwQOASzGDB4A8wAweAJBAg3cQ80WDLAyyMMjCHho8\nALgUM3gAyAPM4AEACTR4BzFfNMjCIAuDLOyhwQOASzGDB4A8wAweAJBAg3cQ80WDLAyyMMjCHho8\nALgUM3gAyAPM4AEACTR4BzFfNMjCIAuDLOyhwQOASzGDB4A8wAweAJBAg3cQ80WDLAyyMMjCHho8\nALgUM3gAyAPM4AEACTR4BzFfNMjCIAuDLOyhwQOASzGDB4A8wAweAJCQtMGPjo7qqaee0uOPP64f\n/ehH+t3vfjfjed3d3dq3b5+am5t1+fLljBTqRswXDbIwyMIgC3uKkh4sKtIjjzyi8vJyXb9+XS0t\nLVq/fr3KysoS5/T392t4eFgnTpzQpUuXdOrUKR06dCjjhQMAkkt6Be/z+VReXi5JKi4u1qpVq/T+\n++9POaevr0+bN2+WJFVWVmpiYkLRaDQz1bpMbW2t0yXkDLIwyMIgC3tSnsGPjIzo7bffVmVl5ZSP\nRyIRBQKBxDoQCCgSiSxchQCAeUmpwV+/fl3Hjh3Trl27VFxcPO14uhtxbp2rBYPBRbu++edcqcfJ\n9f9m4nQ9Tq67urpyqh4n111dXTlVj9PrdM25TXJyclJtbW360pe+pG984xvTjj///PNat26dampq\nJElNTU1qbW2Vz+eb8fXYJmkEg0F+BP0vsjDIwiALY8G3SVqWpZMnT2rVqlUzNndJqqqq0muvvSZJ\nCofDKikpmbW5Yyq+cQ2yMMjCIAt7ku6iuXjxos6ePavVq1friSeekCQ99NBDGh0dlSTV1dVp48aN\nCoVC2rt3r4qLi9XQ0JD5qgEAc0ra4O+++269+OKLc77Io48+umAFLSb8+GmQhUEWBlnYw52sAOBS\nPIsGAPIAz6IBACTQ4B1kZ3+r25CFQRYGWdhDgwcAl2IGDwB5gBk8ACCBBu8g5osGWRhkYZCFPTR4\nAHApZvAAkAeYwQMAEmjwDmK+aJCFQRYGWdhDgwcAl2IGDwB5gBk8ACCBBu8g5osGWRhkYZCFPTR4\nAHApZvAAkAeYwQMAEmjwDmK+aJCFQRYGWdhDgwcAl2IGDwB5gBk8ACCBBu8g5osGWRhkYZCFPTR4\nAHApZvAAkAeYwQMAEmjwDmK+aJCFQRYGWdhTNNcJnZ2dCoVCWr58uY4ePTrt+ODgoNrb27VixQpJ\nUnV1terr6xe+UgBAWuZs8Pfdd5++/vWv69lnn531nLVr16qlpWVBC1sMamtrnS4hZ5CFQRYGWdgz\n54hmzZo1KikpSXpOlt+nBQCkwPYM3uPxKBwOq7m5WYcPH9aVK1cWoq5FgfmiQRYGWRhkYY/tBl9R\nUaGuri4dOXJE27Zt05EjR+b8nFv/0YLBIGvWrG9Znz9/PqfqcXJ9/vz5nKrH6XW6UtoHPzIyora2\nthnfZP1fP/jBD9TW1qZly5bNeJx98ACQPkf2wUej0cQMfmhoSJJmbe4AgOyZcxfN8ePHdeHCBY2P\nj6uhoUE7duxQLBaTJNXV1encuXPq7e1VQUGBvF6vGhsbM160WwSDQXYJ/BdZGGRhkIU9czb4pqam\npMe3bdumbdu2LVhBAICFwbNoACAP8CwaAEACDd5BdrY/uQ1ZGGRhkIU9NHgAcClm8ACQB5jBAwAS\naPAOYr5okIVBFgZZ2EODBwCXYgYPAHmAGTwAIIEG7yDmiwZZGGRhkIU9NHgAcClm8ACQB5jBAwAS\naPAOYr5okIVBFgZZ2EODBwCXYgYPAHmAGTwAIIEG7yDmiwZZGGRhkIU9NHgAcClm8ACQB5jBAwAS\naPAOYr5okIVBFgZZ2EODBwCXYgYPAHmAGTwAIIEG7yDmiwZZGGRhkIU9NHgAcKk5Z/CdnZ0KhUJa\nvny5jh49OuM53d3dGhgYkNfr1Z49e1RRUTHr6zGDB4D0ZWQGf9999+ngwYNJv+jw8LBOnDihxx57\nTKdOnUqrAABAZszZ4NesWaOSkpJZj/f19Wnz5s2SpMrKSk1MTCgajS5chS7GfNEgC4MsDLKwx/YM\nPhKJKBAIJNaBQECRSMTuywIAbCpaiBdJdyv91lOhhfiyOafUW6j9FeOSpNraWknmCmSmdW1tbdLj\nrBfv+qZcqcep9c2P5Uo9Tq/TldKNTiMjI2pra5vxTdbnn39e69atU01NjSSpqalJra2t8vl8M76W\nm99krX/hTX14I5by+aXeQr20c30GKwLgFvN5k9X2FXxVVZXOnDmjmpoahcNhlZSUzNrc3S7dZu3W\nn2Tm49artMWOLAyysGfOBn/8+HFduHBB4+Pjamho0I4dOxSLfXKVWldXp40bNyoUCmnv3r0qLi5W\nQ0NDxot2i+ICa8Ymz5U9gIXAs2hy0NZTIb26e4PTZQDIITyLBgCQsCC7aDA/s80XS72Fi250w6zV\nIAuDLOyhweeg2Zo4b8oCSAcjGgdxZWKQhUEWBlnYQ4MHAJeiwTuI52wYZGGQhUEW9jCDzyOL8c1X\nAPPHPngXYN884H7sgwcAJNDgHcR80SALgywMsrCHBg8ALkWDdxB7fA2yMMjCIAt7aPAA4FI0eAcx\nXzTIwiALgyzsYR+8C7A/HsBM2AfvYuyPB9yDffAAgAQavIOYLxpkYZCFQRb20OABwKVo8A5ij69B\nFgZZGGRhDw0eAFyKBu8g5osGWRhkYZCFPTR4AHApGryDmC8aZGGQhUEW9nAnq4txhyuwuNHgHRQM\nBjN6hTJbE5+p6Tst01nkE7IwyMIeRjQA4FI0eAdxZWKQhUEWBlnYM+eIZmBgQD09PYrH47r//vv1\n4IMPTjk+ODio9vZ2rVixQpJUXV2t+vr6zFQLAEhZ0iv4eDyu06dP6+DBgzp27Jhef/11XblyZdp5\na9euVXt7u9rb22nuaWCPr0EWBlkYZGFP0gY/NDSklStX6o477lBRUZFqamrU19c37bwsP3EYAJCC\npA0+EokoEAgk1n6/X5FIZMo5Ho9H4XBYzc3NOnz48IxX+JgZ80WDLAyyMMjCHtvbJCsqKtTV1SWv\n16tQKKQjR47omWeeWYjaAAA2JL2C9/v9GhsbS6zHxsbk9/unnLNkyRJ5vV5J0oYNGzQ5OamrV68m\n/aK3ztWCweCiXd/8c7a//q1yJY//zcTpepxcd3V15VQ9Tq67urpyqh6n1+lK+iv7YrGYmpqa9LOf\n/Ux+v18HDhxQY2OjysrKEudEo1Hdfvvt8ng8Ghoa0tNPP62Ojo5ZvyC/ss8IBp25iSMXf5WfU1nk\nIrIwyMKYz6/sSzqiKSws1Pe+9z0dOnQosU2yrKxMvb29kqS6ujqdO3dOvb29KigokNfrVWNj4/z/\nBosM37gGWRhkYZCFPfzS7UWo/oU39eGN2LSP84waIHct+BU8MsupHz9z8Rk1/ChukIVBFvbwqAIA\ncCkavIO4MjHIwiALgyzsocEDgEvR4B1kZ3+r25CFQRYGWdhDgwcAl6LBO4j5okEWBlkYZGEPDR4A\nXIoG7yDmiwZZGGRhkIU93OiEhFJv4Yw3O3GHK5CfeFQB5pSLDycDFpv5PKqAEQ0AuBQN3kHMFw2y\nMMjCIAt7aPAA4FI0eAexx9cgC4MsDLKwhwYPAC5Fg3cQ80WDLAyyMMjCHvbBY07sjwfyE/vgMW/s\njweyh33wAIAEGryDmC8aZGGQhUEW9tDgAcClmMFj3upfeFMf3ojNeIw3YIGFNZ8ZPLtoMG/JGvhM\nu24AZBcN3kHBYNC1d+qlu7XSzVmkiywMsrCHBo+MmO3qvv6FN2e5ui+R/o+r/k+QhbE4s1io7cc0\neActxisT5vJA9rCLBgBcas4GPzAwoKamJu3bt0+//e1vZzynu7tb+/btU3Nzsy5fvrzgRboVe3wN\nsjDIwiALe5I2+Hg8rtOnT+vgwYM6duyYXn/9dV25cmXKOf39/RoeHtaJEyf02GOP6dSpUxktGACQ\nmqQNfmhoSCtXrtQdd9yhoqIi1dTUqK+vb8o5fX192rx5sySpsrJSExMTikajmavYRRbjDH42ZGGQ\nhUEW9iRt8JFIRIFAILH2+/2KRCJJzwkEAtPOAQBk34K8yZrlm2Fdg/miQRYGWRhkYU/SbZJ+v19j\nY2OJ9djYmPx+f9rn3Mrn86m/v3++9brK0qVLyeK/yMIgC4MsDJ/Pl/bnJG3wd955p/79739rZGRE\nfr9fb7zxhhobG6ecU1VVpTNnzqimpkbhcFglJSVJC9m0aVPaRQIA0jfnw8ZCoZB6enoUj8d1//33\na/v27ert7ZUk1dXVSZJOnz6tgYEBFRcXq6GhQZ///OczXzkAIKmsP00SAJAd3MkKAC5FgwcAl8rY\nw8YGBgamzO4ffPDBaed0d3drYGBAXq9Xe/bsUUVFRabKcdRcWZw9e1avvPKKLMvSkiVLtHv3bn3u\nc59zqNrMSeV7QvrkBruf/vSn2r9/v6qrq7NcZXakksXg4KB+/etfKxaLqbS0VK2trdkvNAvmymJ8\nfFy/+tWvFI1GFY/H9c1vflNf+9rXnCk2wzo7OxUKhbR8+XIdPXp0xnPS6ptWBsRiMeuHP/yhNTw8\nbH388cfWj3/8Y+tf//rXlHP+9re/Wb/4xS8sy7KscDhsHTx4MBOlOC6VLC5evGhNTExYlmVZoVDI\nlVmkksPN81pbW63Dhw9bf/nLXxyoNPNSyeLq1avW/v37rdHRUcuyLOuDDz5wotSMSyWLF1980frN\nb35jWdYnOXz3u9+1JicnnSg34/7xj39Yb731lvX444/PeDzdvpmREQ2PODBSyeILX/iCli5dKkm6\n6667ptxX4Bap5CBJv//973Xvvfdq+fLlDlSZHalkEQwGVV1dnbhL3K15pJLFpz71Kf3nP/+RJF27\ndk2lpaUqLCx0otyMW7NmjUpKSmY9nm7fzEiD5xEHRipZ3OqPf/yjNmxYmIf955JUvyf6+vq0detW\nSZLH48lqjdmSShbvvfeerl69qqeeeko/+clP9Nprr2W7zKxIJYstW7boypUr+v73v6/m5mbt2rUr\ny1XmjnT7pqNvslrs0Jzi73//u/70pz/pO9/5jtOlOKKnp0cPP/ywPB6PLMta1N8fsVhMly9f1oED\nB/Tkk0/qpZde0nvvved0WY54+eWXVV5erueee07t7e06ffq0rl275nRZjknn/0VG3mTNxCMO8lWq\nf8933nlHzz33nJ588kktW7YsmyVmRSo5vPXWWzp+/Lgk6cMPP9TAwICKiopUVVWV1VozLZUsAoGA\nSktLddttt+m2227TmjVr9M477+gzn/lMtsvNqFSyCIfD2r59uyQlxjnvvvuu7rzzzqzWmgvS7ZsZ\nuYK/9REHk5OTeuONN6b9J62qqkr82JnKIw7yVSpZjI6O6pe//KX27t2rlStXOlRpZqWSw7PPPquO\njg51dHTo3nvv1e7du13X3KXUsrjnnnt08eJFxeNx3bhxQ5cuXVJZWZlDFWdOKll89rOf1fnz5yVJ\n0WhU7777rlasWOFEuY5Lt29m7E5WHnFgzJXFyZMn9de//lWf/vSnJUmFhYU6fPiwkyVnRCrfEzd1\ndnZq06ZNrt0mmUoWr7zyiv785z/L4/Foy5YteuCBB5wsOWPmymJ8fFydnZ0aGxtTPB7X9u3bXfuc\n+OPHj+vChQsaHx+Xz+fTjh07FIvFJM2vb/KoAgBwKe5kBQCXosEDgEvR4AHApWjwAOBSNHgAcCka\nPAC4FA0eAFyKBg8ALvX/7pERkF6cAs0AAAAASUVORK5CYII=\n", "text": [ "" ] } ], "prompt_number": 3 }, { "cell_type": "code", "collapsed": false, "input": [], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 3 } ], "metadata": {} } ] }