|
| 1 | +{ |
| 2 | + "cells": [ |
| 3 | + { |
| 4 | + "cell_type": "code", |
| 5 | + "execution_count": null, |
| 6 | + "metadata": {}, |
| 7 | + "outputs": [], |
| 8 | + "source": [ |
| 9 | + "#\n", |
| 10 | + "# Program 10.2: Brownian motion (brownian.ipynb)\n", |
| 11 | + "# J Wang, Computational modeling and visualization with Python\n", |
| 12 | + "#\n", |
| 13 | + "\n", |
| 14 | + "import numpy as np, random as rnd\n", |
| 15 | + "import matplotlib.pyplot as plt\n", |
| 16 | + "import matplotlib.animation as am\n", |
| 17 | + "%matplotlib notebook\n", |
| 18 | + "\n", |
| 19 | + "class BrownianMotion: # Brownian motion class\n", |
| 20 | + " def __init__(self, N=400, F=0.005, b=0.1, tau=2., h = 0.5):\n", |
| 21 | + " self.N, self.F, self.b, self.tau, self.h = N, F, b, tau, h\n", |
| 22 | + " self.r, self.v = np.zeros((N, 2)), np.zeros((N, 2))\n", |
| 23 | + " self.t = -np.log(np.random.rand(N))*tau # initial kick times\n", |
| 24 | + "\n", |
| 25 | + " def move(self, r, v, dt): # move between kicks\n", |
| 26 | + " edt = np.exp(-self.b*dt)\n", |
| 27 | + " return r + v*(1-edt)/self.b, v*edt\n", |
| 28 | + "\n", |
| 29 | + " def iterate(self): # advance one step\n", |
| 30 | + " r, v, t, h, F = self.r, self.v, self.t, self.h, self.F # alias\n", |
| 31 | + " for i in range(self.N):\n", |
| 32 | + " if t[i] > h: # no kick within current step\n", |
| 33 | + " dt = h # dt= time to end of step\n", |
| 34 | + " else: # suffers kicks before end of step\n", |
| 35 | + " tot, dt = 0., t[i] # tot=time to last kick\n", |
| 36 | + " while t[i] <= h:\n", |
| 37 | + " r[i], v[i] = self.move(r[i], v[i], dt)\n", |
| 38 | + " phi = rnd.random()*2*np.pi # apply kick\n", |
| 39 | + " v[i] += [F*np.cos(phi), F*np.sin(phi)] \n", |
| 40 | + " tot += dt\n", |
| 41 | + " dt = -np.log(rnd.random())*self.tau # sample dt\n", |
| 42 | + " t[i] += dt # next kick\n", |
| 43 | + " dt = h - tot # dt= to end of current step\n", |
| 44 | + " r[i], v[i] = self.move(r[i], v[i], dt) # no kick, just move\n", |
| 45 | + " t[i] -= h\n", |
| 46 | + " \n", |
| 47 | + "def updatefig(*args): # update figure data\n", |
| 48 | + " bm.iterate()\n", |
| 49 | + " plot.set_data(bm.r[:,0], bm.r[:,1]) # update data\n", |
| 50 | + " return [plot] # return plot object\n", |
| 51 | + " \n", |
| 52 | + "bm = BrownianMotion() # create Brownian model\n", |
| 53 | + "fig = plt.figure()\n", |
| 54 | + "plt.subplot(111, aspect='equal')\n", |
| 55 | + "plot = plt.plot(bm.r[:,0], bm.r[:,1], 'o')[0] # create plot object\n", |
| 56 | + "ani = am.FuncAnimation(fig, updatefig, interval=10, blit=True) # animate \n", |
| 57 | + "plt.xlim(-1., 1.), plt.ylim(-1., 1.), plt.show()" |
| 58 | + ] |
| 59 | + }, |
| 60 | + { |
| 61 | + "cell_type": "code", |
| 62 | + "execution_count": null, |
| 63 | + "metadata": { |
| 64 | + "collapsed": true |
| 65 | + }, |
| 66 | + "outputs": [], |
| 67 | + "source": [] |
| 68 | + } |
| 69 | + ], |
| 70 | + "metadata": { |
| 71 | + "kernelspec": { |
| 72 | + "display_name": "Python 3", |
| 73 | + "language": "python", |
| 74 | + "name": "python3" |
| 75 | + }, |
| 76 | + "language_info": { |
| 77 | + "codemirror_mode": { |
| 78 | + "name": "ipython", |
| 79 | + "version": 3 |
| 80 | + }, |
| 81 | + "file_extension": ".py", |
| 82 | + "mimetype": "text/x-python", |
| 83 | + "name": "python", |
| 84 | + "nbconvert_exporter": "python", |
| 85 | + "pygments_lexer": "ipython3", |
| 86 | + "version": "3.6.1" |
| 87 | + } |
| 88 | + }, |
| 89 | + "nbformat": 4, |
| 90 | + "nbformat_minor": 2 |
| 91 | +} |
0 commit comments