Flocking and Emergent Behavior.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The two examples we did before are examples of physical systems with behavior that can be calculated more or less exactly, according to well-defined equations. But what about chaotic or complex systems where we don't necessarily have a precise mathematical description of how the system evolves?\n",
    "\n",
    "Specifically, let's consider systems of large numbers of autonomous agents who move independently but can influence each other. This kind of situation is extremely common in biology and sociology when examining the behavior of animals and humans in groups.\n",
    "\n",
    "<img width=\"500px\" src=\"images/bird_flock.gif\">\n",
    "<img width=\"500px\" src=\"images/traffic.gif\">\n",
    "<img width=\"500px\" src=\"images/fish.gif\">\n",
    "\n",
    "So how can we model this? It would be incredibly difficult to invidually simulate the motivations of each agent (plus, it would require algorithmic knowledge beyond the scope of this two-hour course). We want a simple system that can handle a lot of agents, but still captures the essential qualities of the system. It turns out that this isn't as hard as it sounds! We can abstract away a lot of the complexity behind what drives real organisms and still obtain the same emergent properties when we lump a lot of them together.\n",
    "\n",
    "The simulation we're going to build is a very general model of flocking movement: it could be fish, birds, people who like to move in groups, etc. How do we represent an agent? \n",
    "Let's make a few assumptions. We assume the agents are identical and stateless - they all react uniformly and don't remember anything about their environment. This means we don't have to store any information about the agent itself beyond universal parameters for all agents. This helps us a lot!\n",
    "\n",
    "Now we can simplify the situation down to something very similar to the gravity simulation! Each agent only has to interact with the other agents and environment around it, so we only need to store a velocity and position for each agent.\n",
    "\n",
    "Let \n",
    "$$\\vec{r_i} = \\vec{r_i}(t)$$\n",
    "$$\\vec{v_i} = \\vec{v_i}(t)$$ \n",
    "designate the position and velocity respectively of bird $i$ as a function of time. Just as with the gravity section, we can use a force to produce an acceleration that changes the velocity: the velocity in turn changes the position. We can still use Euler integration to approximate the continuous position and acceleration functions. The theoretical treatment is quite simple!\n",
    "\n",
    "Let's make a representation for a group of birds as a list of such [$\\vec{r}$, $\\vec{v}$] pairs\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "#This list represents two birds: one with position (10,20) and velocity (3,4)\n",
    "#the other has position (50,50) and velocity (-1,-2)\n",
    "birds = [[[10,20],[3,4]],\n",
    "         [[50,50],[-1,-2]]]\n",
    "#This is another two birds\n",
    "morebirds = [[[40,50],[-2,1]],\n",
    "             [[75,25],[5,6]]]\n",
    "\n",
    "#Now we make a list of lists of both sets of birds. This creates two \"flocks\" of birds with two birds in each.\n",
    "#This is the canonical representation which we will work with\n",
    "flocks = [birds, morebirds]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "But wait! Where do these forces come from? Obviously we don't have a direct force acting on birds to flock together. We choose the forces as we develop the model to approximate real world behavior. The difficulty is choosing forces that produce an accurate reflection of reality.\n",
    "\n",
    "Let's start with a really simple model. We have a square box 100x100 units and we put a bunch of birds in it. For now, the birds don't do anything special. They start with a random position and velocity.\n",
    "\n",
    "To start, let's try to keep the birds inside the box! How do we do this? Well, if a bird gets close to a wall it will see the wall and try to turn away so it doesn't crash into the wall. We can simulate this by adding a force that turns the bird away from the wall if it gets too close. Run the boxes below to setup a function that calculate this force for a set of birds.\n",
    "\n",
    "\n",
    "**Be sure to run the next three cells below (by pressing Ctrl-Enter in each cell) or else everything below this won't work!**\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "#Import libraries to help us with animation and calculation\n",
    "import numpy as np\n",
    "from matplotlib import pyplot as plt\n",
    "from matplotlib import animation\n",
    "import random\n",
    "\n",
    "\"\"\"\n",
    "Generates forces that cause birds to avoid running into the walls\n",
    "\"\"\"\n",
    "def avoid_wall_force(flocks):\n",
    "    forces = np.zeros((len(flocks), len(flocks[0]), 2))\n",
    "    #These constants determine how quickly the bird turns as it gets close to the wall\n",
    "    c = 1000.0\n",
    "    d = 50.0\n",
    "    #Determine how far away the bird is when it reacts. We need to set a minimum distance so the bird has time to turn!\n",
    "    minv = max(2, min(sight_radius, 5))\n",
    "\n",
    "    #Inner function to calculate the force on the bird from the walls\n",
    "    #Python lets you nest functions!\n",
    "    #This will calculate the force on a single bird from the walls\n",
    "    def calculate_wall_force(bird):\n",
    "        force = [0,0]\n",
    "        #If the bird is too close to the left edge\n",
    "        if bird[0][0] < minv:\n",
    "            force[0] += c/bird[0][0]**2\n",
    "        #If the bird is too close to the right edge\n",
    "        if 100 - bird[0][0] < minv:\n",
    "            force[0] -= c/(100-bird[0][0])**2\n",
    "        #If the bird is too close to the bottom edge\n",
    "        if bird[0][1] < minv:\n",
    "            force[1] += c/bird[0][1]**2\n",
    "        #If the bird is too close to the top edge\n",
    "        if 100 - bird[0][1] < minv:\n",
    "            force[1] -= c/(100-bird[0][1])**2\n",
    "        return force\n",
    "\n",
    "    #Update the force matrix with the forces we calculated above\n",
    "    for i in xrange(len(flocks)):\n",
    "        forces[i] = map(calculate_wall_force, flocks[i])\n",
    "    \n",
    "    return forces"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Birds also like to keep moving, so let's also add a force that will try to keep them moving at around a constant fraction of their max speed if no other forces are acting on them."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "\"\"\"Birds always want to be in motion: this force will drive them to try to maintain a constant flight speed of p*max_vel\"\"\"\n",
    "def motive_force(flocks):\n",
    "    forces = np.zeros((len(flocks), len(flocks[0]), 2))\n",
    "    c = 10.0\n",
    "    p = 0.75\n",
    "    for i in xrange(len(flocks)):\n",
    "        forces[i] = map(lambda bird: np.multiply(c*p*bird_max_vel/np.linalg.norm(bird[1]), bird[1]), flocks[i])\n",
    "\n",
    "    return forces"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Let's also put some boilerplate code up here to help us generate and render the birds correctly. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "#Visibility flag to help us render the birds\n",
    "are_visible = 0\n",
    "\n",
    "\"\"\"This calculates the sum of an arbitrary number of arrays. Used for totalling the force arrays later\"\"\"\n",
    "def sum_arrays(*arrays):\n",
    "    if len(arrays) == 1:\n",
    "        return arrays[0]\n",
    "    total = arrays[0]\n",
    "    for array in arrays[1:]:\n",
    "        total = np.add(total, array)\n",
    "    return total\n",
    "\n",
    "\"\"\"This function updates the velocity and position of a set of agents based on the forces that are applied\n",
    "'patches' are just rendering objects used by matplotlib: don't worry about what happens to them\"\"\"\n",
    "def updatePositionVelocity(agents, forces, patches, npatches):\n",
    "    for i,flock in enumerate(agents):\n",
    "        #Update position and velocities using velocity and force: assume birds have mass 1\n",
    "        positions = np.add(map(lambda bird: bird[0], flock), np.multiply(interval/1000.0,map(lambda bird: bird[1], flock)))\n",
    "        velocities = np.add(map(lambda bird: bird[1], flock), np.multiply(interval/1000.0,forces[i]))\n",
    "        #Cap velocities at max_vel\n",
    "        velocities = map(lambda velocity: velocity if np.linalg.norm(velocity) < bird_max_vel else np.divide(velocity, np.linalg.norm(velocity)/bird_max_vel), velocities)\n",
    "        #Update our flocks with the new positions\n",
    "        agents[i] = map(lambda position, velocity: (position, velocity), positions, velocities)\n",
    "        newpolys = map(lambda bird: gen_polygon(bird[0], bird[1]), flock)\n",
    "        for j,poly in enumerate(newpolys):\n",
    "            patches[i][j].set_xy(poly)\n",
    "            npatches.append(patches[i][j])\n",
    "\n",
    "    return npatches\n",
    "\n",
    "\"\"\" This function will randomly initialize the birds. \n",
    "It produces num_flocks different flocks each with num_bird_in_flock birds\"\"\"\n",
    "def random_init_birds(num_flocks, num_birds_in_flock):\n",
    "    if num_flocks > len(flock_colors):\n",
    "        raise \"Too many flocks!, Add more colors!\"\n",
    "    flocks = []\n",
    "    for i in xrange(num_flocks):\n",
    "        flocks.append([])\n",
    "        for j in xrange(num_birds_in_flock):\n",
    "            #Each bird is represented by a pair of vectors [position, velocity]\n",
    "            #We randomly generate these positions and vectors\n",
    "            flocks[i].append([(95*random.random()+5, 95*random.random()+5), generate_random_velocity_vector(bird_max_vel)])\n",
    "    return flocks\n",
    "\n",
    "#Generates a random velocity vector with magnitude at most max_vel\n",
    "def generate_random_velocity_vector(max_vel):\n",
    "    x,y = max_vel, max_vel\n",
    "    while np.linalg.norm([x,y]) > max_vel:\n",
    "        x = (random.random()-0.5)*2*max_vel\n",
    "        y = (random.random()-0.5)*2*max_vel\n",
    "    return (x,y)\n",
    "\n",
    "\n",
    "\"\"\"Generates the icon for a single bird on the plot. Don't worry about this stuff ;)\"\"\"\n",
    "def gen_polygon(x, v):\n",
    "    #Scales the size of a bird icon\n",
    "    c = 2.0\n",
    "    if x[0] > 100 or x[0] < 0 or x[1] > 100 or x[1] < 0:\n",
    "        raise ValueError(\"Bird out of bounds at \" + str(x[0]) + \",\" + str(x[1]))\n",
    "    vn = np.linalg.norm(v)\n",
    "    return [[x[0] + (v[0]/vn)*(2*c/3), x[1] + (v[1]/vn)*(2*c/3)],\n",
    "            [x[0] - (v[0]/vn)*(c/3) + (v[1]/vn)*(c/4),   x[1] - (v[1]/vn)*(2*c/3) - (v[0]/vn)*(c/4)],\n",
    "            [x[0] - (v[0]/vn)*(c/3) - (v[1]/vn)*(c/4),   x[1] - (v[1]/vn)*(2*c/3) + (v[0]/vn)*(c/4)]]\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Alright, now that we've set up the infrastructure, let's go ahead and run our first simulation!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "#Setup figure\n",
    "#Don't change this!\n",
    "fig = plt.figure()\n",
    "ax = plt.axes(xlim=(0, 100), ylim=(0, 100))\n",
    "plt.tight_layout()\n",
    "\n",
    "#Colors\n",
    "#Feel free to change!\n",
    "flock_colors = ['b','g','c','grey', 'darkviolet','gold']\n",
    "\n",
    "#Set number of blocks and number of birds in each flock\n",
    "#Don't set these too high or else perfomance will suffer!\n",
    "#These need to ints\n",
    "num_flocks = 1\n",
    "num_birds_in_flock = 20\n",
    "\n",
    "#Flocking settings -try changing these and see what happens! \n",
    "#ALL OF THESE SHOULD BE FLOATS!\n",
    "bird_max_vel = 25.0\n",
    "sight_radius = 10.0\n",
    "\n",
    "#Rendering speeds\n",
    "interval = 20.0\n",
    "\n",
    "#Populates and stores birds in the simulation\n",
    "#Each flock is a row of birds, each bird is a pair of vectors [position, velocity]\n",
    "flocks = random_init_birds(num_flocks, num_birds_in_flock)\n",
    "patches = [[plt.Polygon(gen_polygon(flocks[i][j][0], flocks[i][j][1]), alpha=1, color=flock_colors[i]) for j in xrange(len(flocks[i]))] for i in xrange(len(flocks))]\n",
    "\n",
    "def init():\n",
    "    return ax\n",
    "\n",
    "#Update the frame\n",
    "def animate(i):\n",
    "    npatches = []\n",
    "    global flocks\n",
    "    global patches\n",
    "\n",
    "    #If this is the first frame, initialize the view!\n",
    "    if i==0:\n",
    "        #Add all of our birds to the window\n",
    "        for flock in patches:\n",
    "            for patch in flock:\n",
    "                plt.gca().add_patch(patch)\n",
    "                npatches.append(patch)\n",
    "        return npatches\n",
    "\n",
    "    #Calculate the wall avoidance force\n",
    "    wall_force = avoid_wall_force(flocks)\n",
    "    #Calculate the motive force that tries to keep the birds moving\n",
    "    mot_force = motive_force(flocks)\n",
    "    #Sum the forces to obtain the total force on each bird!\n",
    "    total_force = sum_arrays(wall_force,mot_force)\n",
    "\n",
    "    #Now we update the position and velocity of the birds according to the force\n",
    "    updatePositionVelocity(flocks, total_force, patches, npatches)\n",
    "\n",
    "    return npatches\n",
    "\n",
    "#Animate!\n",
    "anim = animation.FuncAnimation(fig, animate, \n",
    "                               init_func=init, \n",
    "                               frames=10000, \n",
    "                               interval=interval,\n",
    "                               blit=True)\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Well that wasn't very interesting. The birds just ignored each other and went on more or less straight paths, but at least we managed to keep them in the box! Now we can add some more behaviors. Let's think about the kinds of behaviors that define a flock. What simple forces and rules can we add to produce an emergent system that acts like a flock?\n",
    "\n",
    "It turns out there is a subfield of scientists that research just these kinds of concerns, and flocking is extremely well studied.\n",
    "\n",
    "There are only three individual \"fundamental forces\" that combine produce flocking behavior!\n",
    "\n",
    "For realism, let's introduce the idea of a 'sight range' and 'separation distance.' Birds can only directly react to things within their sight range, and will try to stay away from other birds within the separation distance. The separation distance also helps the visualization, since it's easier to tell what's going on when all the birds aren't overlapping.\n",
    "\n",
    "<img width=\"225px\" src=\"images/detection_circles.gif\">\n",
    "\n",
    "The three fundamental flocking forces are [drumroll please] . . . \n",
    "\n",
    "### Cohesion\n",
    "\n",
    "In order for a flock to be, well, a flock, the birds need to stay close together! The cohesion force promotes staying together by adding a force that pushes a bird toward the average location of the birds within its *sight range* (i.e. the center of the visible flock). To produce this force, we simply average the positions of the birds that a single bird can see and make a force vector that points toward that average.\n",
    "\n",
    "<img width=\"225px\" src=\"images/cohesion.gif\">\n",
    "\n",
    "### Alignment\n",
    "\n",
    "Flocks move in a coordinated manner, as is clearly visible from the gifs in the introduction. This coordination is produced by every member of the flock trying to stay in alignment with its neighbors. We can simulate this by adding a force that pushes the velocity of a bird towards the average velocity of all the birds around it.\n",
    "\n",
    "<img width=\"225px\" src=\"images/alignment.gif\">\n",
    "\n",
    "### Separation\n",
    "\n",
    "Wait, didn't we just say flocks needed to say together? Well, yes, but of course birds (and flocking animals in general) don't want to be right on top of each other. We all need personal space! We can simulate this by adding a force that pushes a bird away from its neighbors that are within the separation distance.\n",
    "\n",
    "<img width=\"225px\" src=\"images/separation.gif\">\n",
    "\n",
    "\n",
    "\n",
    "We now have a set of forces that can produce flocking. But wait, you might say! While we know the general form of these forces, we don't have exact descriptions! Should the force scale with distance, or the square of distance, or not at all? And how strong should the alignment force be compared to the cohesion and separation forces?\n",
    "\n",
    "The answer is we don't precisely know! We tune the parameters to make the model behave realistically. One possible setting of these parameters is used below, but we encourage you to play around with the parameters in the section below and see the results!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**The next box defines functions to calculate the three forces we described above. Run it before proceeding.**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "\"\"\"\n",
    "Generates forces that tend to bring birds closer together, at least until they get too close for comfort\n",
    "\"\"\"\n",
    "def cohesion_force(flocks):\n",
    "    forces = np.zeros((len(flocks), len(flocks[0]), 2))\n",
    "    c = 1.0\n",
    "\n",
    "    for i in xrange(len(flocks)):\n",
    "        counts = [1 for k in xrange(len(flocks[i]))]\n",
    "        for j in xrange(len(flocks[i])):\n",
    "            for k in xrange(j+1, len(flocks[i])):\n",
    "                #Calculate offset vector between the two birds\n",
    "                f = np.subtract(flocks[i][j][0], flocks[i][k][0])\n",
    "                dist = np.linalg.norm(f)\n",
    "                #If the birds are close enough to see each other but not too close, push them towards each other\n",
    "                if dist > individual_separation and dist < sight_radius:\n",
    "                    forces[i][j] = np.add(forces[i][j], np.multiply(-cohesion*c, f))\n",
    "                    forces[i][k] = np.add(forces[i][k], np.multiply(cohesion*c, f))\n",
    "                    counts[i] += 1\n",
    "                    counts[k] += 1\n",
    "        forces[i] = np.divide(forces[i],np.transpose([counts, counts]))\n",
    "\n",
    "    return forces\n",
    "\n",
    "\"\"\"\n",
    "Generates a force that aligns the movement of nearby members of the flock\n",
    "\"\"\"\n",
    "def alignment_force(flocks):\n",
    "    forces = np.zeros((len(flocks), len(flocks[0]), 2))\n",
    "    c = 10.0\n",
    "\n",
    "    for i in xrange(len(flocks)):\n",
    "        for j in xrange(len(flocks[i])):\n",
    "            for k in xrange(j+1, len(flocks[i])):\n",
    "                dist = np.linalg.norm(np.subtract(flocks[i][j][0], flocks[i][k][0]))\n",
    "                #If the birds are close enough to see each other, push them towards alignment\n",
    "                if dist < sight_radius:\n",
    "                    #Force vector is proportional to the difference in velocity vectors\n",
    "                    f = np.subtract(flocks[i][k][1], flocks[i][j][1])\n",
    "                    forces[i][j] = np.add(forces[i][j], np.multiply(alignment*c/dist, f))\n",
    "                    forces[i][k] = np.add(forces[i][k], np.multiply(-alignment*c/dist, f))\n",
    "\n",
    "    return forces\n",
    "\n",
    "\"\"\" Birds don't want to be too close together! The separation force drives them apart if they are closer\n",
    "than individual_separation apart\"\"\"\n",
    "def separation_force(flocks):\n",
    "    forces = np.zeros((len(flocks), len(flocks[0]), 2))\n",
    "    c =500.0\n",
    "\n",
    "    for i in xrange(len(flocks)):\n",
    "        for j in xrange(len(flocks[i])):\n",
    "            for k in xrange(j+1, len(flocks[i])):\n",
    "                f = np.subtract(flocks[i][k][0], flocks[i][j][0])\n",
    "                dist = np.linalg.norm(f)\n",
    "                #If the birds are too close, push them apart\n",
    "                if dist < individual_separation:\n",
    "                    forces[i][j] = np.add(forces[i][j], np.multiply(-c/dist, f))\n",
    "                    forces[i][k] = np.add(forces[i][k], np.multiply(c/dist, f))\n",
    "\n",
    "    return forces"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now let's do that simulation again, but with the extra forces we just added.\n",
    "\n",
    "Note that only birds of the same color will flock with each other! (These birds are kind of discriminatory).\n",
    "\n",
    "Play with the parameters defined in the top section of the code and see how they affect the flocks that form. In particular, check out how the flock size and individual_separation interact as the flocks get larger."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "ename": "ValueError",
     "evalue": "Bird out of bounds at 58.1544811856,100.024081197",
     "output_type": "error",
     "traceback": [
      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[1;31mValueError\u001b[0m                                Traceback (most recent call last)",
      "\u001b[1;32mC:\\Users\\Lei\\AppData\\Local\\Enthought\\Canopy\\User\\lib\\site-packages\\matplotlib\\backends\\backend_wx.pyc\u001b[0m in \u001b[0;36m_on_timer\u001b[1;34m(self, *args)\u001b[0m\n\u001b[0;32m    210\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    211\u001b[0m     \u001b[1;32mdef\u001b[0m \u001b[0m_on_timer\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mself\u001b[0m\u001b[1;33m,\u001b[0m \u001b[1;33m*\u001b[0m\u001b[0margs\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m--> 212\u001b[1;33m         \u001b[0mTimerBase\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0m_on_timer\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mself\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m    213\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    214\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;32mC:\\Users\\Lei\\AppData\\Local\\Enthought\\Canopy\\User\\lib\\site-packages\\matplotlib\\backend_bases.pyc\u001b[0m in \u001b[0;36m_on_timer\u001b[1;34m(self)\u001b[0m\n\u001b[0;32m   1271\u001b[0m         '''\n\u001b[0;32m   1272\u001b[0m         \u001b[1;32mfor\u001b[0m \u001b[0mfunc\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0margs\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mkwargs\u001b[0m \u001b[1;32min\u001b[0m \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mcallbacks\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m-> 1273\u001b[1;33m             \u001b[0mret\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mfunc\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;33m*\u001b[0m\u001b[0margs\u001b[0m\u001b[1;33m,\u001b[0m \u001b[1;33m**\u001b[0m\u001b[0mkwargs\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m   1274\u001b[0m             \u001b[1;31m# docstring above explains why we use `if ret == False` here,\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m   1275\u001b[0m             \u001b[1;31m# instead of `if not ret`.\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;32mC:\\Users\\Lei\\AppData\\Local\\Enthought\\Canopy\\User\\lib\\site-packages\\matplotlib\\animation.pyc\u001b[0m in \u001b[0;36m_step\u001b[1;34m(self, *args)\u001b[0m\n\u001b[0;32m    908\u001b[0m         \u001b[1;31m# delay and set the callback to one which will then set the interval\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    909\u001b[0m         \u001b[1;31m# back.\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m--> 910\u001b[1;33m         \u001b[0mstill_going\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mAnimation\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0m_step\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mself\u001b[0m\u001b[1;33m,\u001b[0m \u001b[1;33m*\u001b[0m\u001b[0margs\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m    911\u001b[0m         \u001b[1;32mif\u001b[0m \u001b[1;32mnot\u001b[0m \u001b[0mstill_going\u001b[0m \u001b[1;32mand\u001b[0m \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mrepeat\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    912\u001b[0m             \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mframe_seq\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mnew_frame_seq\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;32mC:\\Users\\Lei\\AppData\\Local\\Enthought\\Canopy\\User\\lib\\site-packages\\matplotlib\\animation.pyc\u001b[0m in \u001b[0;36m_step\u001b[1;34m(self, *args)\u001b[0m\n\u001b[0;32m    767\u001b[0m         \u001b[1;32mtry\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    768\u001b[0m             \u001b[0mframedata\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mnext\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mframe_seq\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m--> 769\u001b[1;33m             \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0m_draw_next_frame\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mframedata\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0m_blit\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m    770\u001b[0m             \u001b[1;32mreturn\u001b[0m \u001b[0mTrue\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    771\u001b[0m         \u001b[1;32mexcept\u001b[0m \u001b[0mStopIteration\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;32mC:\\Users\\Lei\\AppData\\Local\\Enthought\\Canopy\\User\\lib\\site-packages\\matplotlib\\animation.pyc\u001b[0m in \u001b[0;36m_draw_next_frame\u001b[1;34m(self, framedata, blit)\u001b[0m\n\u001b[0;32m    786\u001b[0m         \u001b[1;31m# post- draw, as well as the drawing of the frame itself.\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    787\u001b[0m         \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0m_pre_draw\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mframedata\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mblit\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m--> 788\u001b[1;33m         \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0m_draw_frame\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mframedata\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m    789\u001b[0m         \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0m_post_draw\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mframedata\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mblit\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    790\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;32mC:\\Users\\Lei\\AppData\\Local\\Enthought\\Canopy\\User\\lib\\site-packages\\matplotlib\\animation.pyc\u001b[0m in \u001b[0;36m_draw_frame\u001b[1;34m(self, framedata)\u001b[0m\n\u001b[0;32m   1089\u001b[0m         \u001b[1;31m# Call the func with framedata and args. If blitting is desired,\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m   1090\u001b[0m         \u001b[1;31m# func needs to return a sequence of any artists that were modified.\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m-> 1091\u001b[1;33m         \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0m_drawn_artists\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0m_func\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mframedata\u001b[0m\u001b[1;33m,\u001b[0m \u001b[1;33m*\u001b[0m\u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0m_args\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m",
      "\u001b[1;32m<ipython-input-11-05c6358468ed>\u001b[0m in \u001b[0;36manimate\u001b[1;34m(i)\u001b[0m\n\u001b[0;32m     63\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m     64\u001b[0m     \u001b[1;31m#Now we update the position and velocity of the birds according to the force\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m---> 65\u001b[1;33m     \u001b[0mupdatePositionVelocity\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mflocks\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mtotal_force\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mpatches\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mnpatches\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m     66\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m     67\u001b[0m     \u001b[1;32mreturn\u001b[0m \u001b[0mnpatches\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;32m<ipython-input-4-3fc8f8aa9ba2>\u001b[0m in \u001b[0;36mupdatePositionVelocity\u001b[1;34m(agents, forces, patches, npatches)\u001b[0m\n\u001b[0;32m     22\u001b[0m         \u001b[1;31m#Update our flocks with the new positions\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m     23\u001b[0m         \u001b[0magents\u001b[0m\u001b[1;33m[\u001b[0m\u001b[0mi\u001b[0m\u001b[1;33m]\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mmap\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;32mlambda\u001b[0m \u001b[0mposition\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mvelocity\u001b[0m\u001b[1;33m:\u001b[0m \u001b[1;33m(\u001b[0m\u001b[0mposition\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mvelocity\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mpositions\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mvelocities\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m---> 24\u001b[1;33m         \u001b[0mnewpolys\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mmap\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;32mlambda\u001b[0m \u001b[0mbird\u001b[0m\u001b[1;33m:\u001b[0m \u001b[0mgen_polygon\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mbird\u001b[0m\u001b[1;33m[\u001b[0m\u001b[1;36m0\u001b[0m\u001b[1;33m]\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mbird\u001b[0m\u001b[1;33m[\u001b[0m\u001b[1;36m1\u001b[0m\u001b[1;33m]\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mflock\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m     25\u001b[0m         \u001b[1;32mfor\u001b[0m \u001b[0mj\u001b[0m\u001b[1;33m,\u001b[0m\u001b[0mpoly\u001b[0m \u001b[1;32min\u001b[0m \u001b[0menumerate\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mnewpolys\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m     26\u001b[0m             \u001b[0mpatches\u001b[0m\u001b[1;33m[\u001b[0m\u001b[0mi\u001b[0m\u001b[1;33m]\u001b[0m\u001b[1;33m[\u001b[0m\u001b[0mj\u001b[0m\u001b[1;33m]\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mset_xy\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mpoly\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;32m<ipython-input-4-3fc8f8aa9ba2>\u001b[0m in \u001b[0;36m<lambda>\u001b[1;34m(bird)\u001b[0m\n\u001b[0;32m     22\u001b[0m         \u001b[1;31m#Update our flocks with the new positions\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m     23\u001b[0m         \u001b[0magents\u001b[0m\u001b[1;33m[\u001b[0m\u001b[0mi\u001b[0m\u001b[1;33m]\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mmap\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;32mlambda\u001b[0m \u001b[0mposition\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mvelocity\u001b[0m\u001b[1;33m:\u001b[0m \u001b[1;33m(\u001b[0m\u001b[0mposition\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mvelocity\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mpositions\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mvelocities\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m---> 24\u001b[1;33m         \u001b[0mnewpolys\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mmap\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;32mlambda\u001b[0m \u001b[0mbird\u001b[0m\u001b[1;33m:\u001b[0m \u001b[0mgen_polygon\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mbird\u001b[0m\u001b[1;33m[\u001b[0m\u001b[1;36m0\u001b[0m\u001b[1;33m]\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mbird\u001b[0m\u001b[1;33m[\u001b[0m\u001b[1;36m1\u001b[0m\u001b[1;33m]\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mflock\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m     25\u001b[0m         \u001b[1;32mfor\u001b[0m \u001b[0mj\u001b[0m\u001b[1;33m,\u001b[0m\u001b[0mpoly\u001b[0m \u001b[1;32min\u001b[0m \u001b[0menumerate\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mnewpolys\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m     26\u001b[0m             \u001b[0mpatches\u001b[0m\u001b[1;33m[\u001b[0m\u001b[0mi\u001b[0m\u001b[1;33m]\u001b[0m\u001b[1;33m[\u001b[0m\u001b[0mj\u001b[0m\u001b[1;33m]\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mset_xy\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mpoly\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;32m<ipython-input-4-3fc8f8aa9ba2>\u001b[0m in \u001b[0;36mgen_polygon\u001b[1;34m(x, v)\u001b[0m\n\u001b[0;32m     57\u001b[0m     \u001b[0mc\u001b[0m \u001b[1;33m=\u001b[0m \u001b[1;36m2.0\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m     58\u001b[0m     \u001b[1;32mif\u001b[0m \u001b[0mx\u001b[0m\u001b[1;33m[\u001b[0m\u001b[1;36m0\u001b[0m\u001b[1;33m]\u001b[0m \u001b[1;33m>\u001b[0m \u001b[1;36m100\u001b[0m \u001b[1;32mor\u001b[0m \u001b[0mx\u001b[0m\u001b[1;33m[\u001b[0m\u001b[1;36m0\u001b[0m\u001b[1;33m]\u001b[0m \u001b[1;33m<\u001b[0m \u001b[1;36m0\u001b[0m \u001b[1;32mor\u001b[0m \u001b[0mx\u001b[0m\u001b[1;33m[\u001b[0m\u001b[1;36m1\u001b[0m\u001b[1;33m]\u001b[0m \u001b[1;33m>\u001b[0m \u001b[1;36m100\u001b[0m \u001b[1;32mor\u001b[0m \u001b[0mx\u001b[0m\u001b[1;33m[\u001b[0m\u001b[1;36m1\u001b[0m\u001b[1;33m]\u001b[0m \u001b[1;33m<\u001b[0m \u001b[1;36m0\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m---> 59\u001b[1;33m         \u001b[1;32mraise\u001b[0m \u001b[0mValueError\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;34m\"Bird out of bounds at \"\u001b[0m \u001b[1;33m+\u001b[0m \u001b[0mstr\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mx\u001b[0m\u001b[1;33m[\u001b[0m\u001b[1;36m0\u001b[0m\u001b[1;33m]\u001b[0m\u001b[1;33m)\u001b[0m \u001b[1;33m+\u001b[0m \u001b[1;34m\",\"\u001b[0m \u001b[1;33m+\u001b[0m \u001b[0mstr\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mx\u001b[0m\u001b[1;33m[\u001b[0m\u001b[1;36m1\u001b[0m\u001b[1;33m]\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m     60\u001b[0m     \u001b[0mvn\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mnp\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mlinalg\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mnorm\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mv\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m     61\u001b[0m     return [[x[0] + (v[0]/vn)*(2*c/3), x[1] + (v[1]/vn)*(2*c/3)],\n",
      "\u001b[1;31mValueError\u001b[0m: Bird out of bounds at 58.1544811856,100.024081197"
     ]
    }
   ],
   "source": [
    "#Setup figure\n",
    "#Don't change this!\n",
    "fig = plt.figure()\n",
    "ax = plt.axes(xlim=(0, 100), ylim=(0, 100))\n",
    "plt.tight_layout()\n",
    "\n",
    "#Colors\n",
    "#Feel free to change!\n",
    "flock_colors = ['b','g','c','grey', 'darkviolet','gold']\n",
    "\n",
    "#Set number of blocks and number of birds in each flock\n",
    "#Don't set these too high or else perfomance will suffer!\n",
    "#These need to ints\n",
    "num_flocks = 2\n",
    "num_birds_in_flock = 15\n",
    "\n",
    "#Flocking settings -try changing these and see what happens! \n",
    "#ALL OF THESE SHOULD BE FLOATS!\n",
    "bird_max_vel = 25.0\n",
    "sight_radius = 15.0\n",
    "individual_separation = 1.0\n",
    "cohesion = 1.0\n",
    "alignment = 1.0\n",
    "\n",
    "#Rendering speeds\n",
    "interval = 20.0\n",
    "\n",
    "#Populates and stores birds in the simulation\n",
    "#Each flock is a row of birds, each bird is a pair of vectors [position, velocity]\n",
    "flocks = random_init_birds(num_flocks, num_birds_in_flock)\n",
    "patches = [[plt.Polygon(gen_polygon(flocks[i][j][0], flocks[i][j][1]), alpha=1, color=flock_colors[i]) for j in xrange(len(flocks[i]))] for i in xrange(len(flocks))]\n",
    "\n",
    "def init():\n",
    "    return ax\n",
    "\n",
    "#Update the frame\n",
    "def animate(i):\n",
    "    npatches = []\n",
    "    global flocks\n",
    "    global patches\n",
    "\n",
    "    #If this is the first frame, initialize the view!\n",
    "    if i==0:\n",
    "        #Add all of our birds to the window\n",
    "        for flock in patches:\n",
    "            for patch in flock:\n",
    "                plt.gca().add_patch(patch)\n",
    "                npatches.append(patch)\n",
    "        return npatches\n",
    "\n",
    "    #Calculate the wall avoidance force\n",
    "    wall_force = avoid_wall_force(flocks)\n",
    "    #Calculate the motive force that tries to keep the birds moving\n",
    "    mot_force = motive_force(flocks)\n",
    "    #Calculate the cohesion force\n",
    "    coh_force = cohesion_force(flocks)\n",
    "    #Calculate the alignment force\n",
    "    align_force = alignment_force(flocks)\n",
    "    #Calculate the separation force\n",
    "    sep_force = separation_force(flocks)\n",
    "    #Sum the forces to obtain the total force on each bird!\n",
    "    total_force = sum_arrays(wall_force,mot_force, coh_force, align_force, sep_force)\n",
    "\n",
    "    #Now we update the position and velocity of the birds according to the force\n",
    "    updatePositionVelocity(flocks, total_force, patches, npatches)\n",
    "\n",
    "    return npatches\n",
    "\n",
    "#Animate!\n",
    "anim = animation.FuncAnimation(fig, animate, \n",
    "                               init_func=init, \n",
    "                               frames=10000, \n",
    "                               interval=interval,\n",
    "                               blit=True)\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally, just for fun let's add a hawk into the simulation\n",
    "\n",
    "The birds all try to run away from the hawk, while the hawk chases the closest bird that it can see.\n",
    "We simulate this by adding a force that pushes the birds away from the hawk if the hawk is within the bird's sight radius, while adding a force to the hawk that aims it toward the closest bird it can see.\n",
    "\n",
    "If you have time, try changing the hawk AI so it is smarter when pursuing a bird (maybe it tries to predict the bird's position according to its velocity instead of just going towards its current position) or implementing a way for the hawk to \"eat\" [remove from the field] birds that it gets close enough to.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "ename": "NameError",
     "evalue": "name 'plt' is not defined",
     "output_type": "error",
     "traceback": [
      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[1;31mNameError\u001b[0m                                 Traceback (most recent call last)",
      "\u001b[1;32m<ipython-input-1-3b9efd70ec3d>\u001b[0m in \u001b[0;36m<module>\u001b[1;34m()\u001b[0m\n\u001b[0;32m      1\u001b[0m \u001b[1;31m#Setup figure\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m      2\u001b[0m \u001b[1;31m#Don't change this!\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m----> 3\u001b[1;33m \u001b[0mfig\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mplt\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mfigure\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m      4\u001b[0m \u001b[0max\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mplt\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0maxes\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mxlim\u001b[0m\u001b[1;33m=\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;36m0\u001b[0m\u001b[1;33m,\u001b[0m \u001b[1;36m100\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mylim\u001b[0m\u001b[1;33m=\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;36m0\u001b[0m\u001b[1;33m,\u001b[0m \u001b[1;36m100\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m      5\u001b[0m \u001b[0mplt\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mtight_layout\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;31mNameError\u001b[0m: name 'plt' is not defined"
     ]
    }
   ],
   "source": [
    "#Setup figure\n",
    "#Don't change this!\n",
    "fig = plt.figure()\n",
    "ax = plt.axes(xlim=(0, 100), ylim=(0, 100))\n",
    "plt.tight_layout()\n",
    "\n",
    "#Colors\n",
    "#Feel free to change!\n",
    "flock_colors = ['b','g','c','grey', 'darkviolet','gold']\n",
    "pred_color = 'r'\n",
    "\n",
    "#Set number of blocks and number of birds in each flock\n",
    "#Don't set these too high or else perfomance will suffer!\n",
    "#These need to ints\n",
    "num_flocks = 4\n",
    "num_birds_in_flock = 15\n",
    "\n",
    "#Predator settings\n",
    "#Note: With the way the simulation is set up, the chases tend to most interesting when the predator is a little slower\n",
    "#than the birds but can see somewhat further\n",
    "pred_sight_radius = 20.0\n",
    "pred_max_vel = 20.0\n",
    "#Set pred_chase to zero to have the predator ignore the birds\n",
    "pred_chase = 1\n",
    "\n",
    "#Flocking settings -try changing these and see what happens! \n",
    "#ALL OF THESE SHOULD BE FLOATS!\n",
    "bird_max_vel = 25.0\n",
    "sight_radius = 10.0\n",
    "individual_separation = 1.0\n",
    "cohesion = 1.0\n",
    "alignment = 1.0\n",
    "predator_avoidance = 1.0\n",
    "\n",
    "#Rendering speeds\n",
    "interval = 20.0\n",
    "\n",
    "#Populates and stores birds in the simulation\n",
    "#Each flock is a row of birds, each bird is a pair of vectors [position, velocity]\n",
    "flocks = random_init_birds(num_flocks, num_birds_in_flock)\n",
    "patches = [[plt.Polygon(gen_polygon(flocks[i][j][0], flocks[i][j][1]), alpha=1, color=flock_colors[i]) for j in xrange(len(flocks[i]))] for i in xrange(len(flocks))]\n",
    "#Oh no it's a hawk!\n",
    "predator = [[95*random.random()+5, 95*random.random()+5], generate_random_velocity_vector(bird_max_vel)]\n",
    "pred_patch = plt.Polygon(gen_polygon(predator[0], predator[1]), alpha=1, color=pred_color)\n",
    "\n",
    "\n",
    "\"\"\"Produces the effect of a predator on the motion of the birds and on the motion of the predator itself\n",
    "The predator will track towards the closest bird in its sight radius\n",
    "The birds always try to flee the predator\"\"\"\n",
    "def predator_force(flocks, predator):\n",
    "    forces = np.zeros((len(flocks), len(flocks[0]), 2))\n",
    "    pred_force = [0,0]\n",
    "    c = 500.0\n",
    "\n",
    "    min_distance_bird = [100, None]\n",
    "    for i in xrange(len(flocks)):\n",
    "        for j in xrange(len(flocks[i])):\n",
    "            f = np.subtract(predator[0], flocks[i][j][0])\n",
    "            dist = np.linalg.norm(f)\n",
    "            if dist < min_distance_bird[0] and dist < pred_sight_radius:\n",
    "                min_distance_bird = [dist, -f]\n",
    "            if dist < sight_radius:\n",
    "                forces[i][j] = np.add(forces[i][j], np.multiply(-predator_avoidance*c/dist, f))\n",
    "\n",
    "    if min_distance_bird[1] != None and pred_chase > 0:\n",
    "        pred_force = np.multiply(10.0,min_distance_bird[1])\n",
    "\n",
    "    return pred_force,forces\n",
    "\n",
    "#Helper function to update the predator's position\n",
    "def updatePositionVelocityForPredator(predator, force, pred_patch, npatches):\n",
    "    pred_patch.set_xy(gen_polygon(predator[0], predator[1]))\n",
    "    predator[0] = np.add(predator[0], np.multiply(interval/1000.0, predator[1]))\n",
    "    predator[1] = np.add(predator[1], np.multiply(interval/1000.0, force))\n",
    "    #Cap predator velocity also!\n",
    "    predator[1] = predator[1] if np.linalg.norm(predator[1]) < pred_max_vel else np.divide(predator[1], np.linalg.norm(predator[1])/pred_max_vel)\n",
    "    npatches.append(pred_patch)\n",
    "\n",
    "    return npatches\n",
    "\n",
    "\n",
    "def init():\n",
    "    return ax,\n",
    "\n",
    "#Update the frame\n",
    "def animate(i):\n",
    "    npatches = []\n",
    "    global flocks\n",
    "    global patches\n",
    "    global predator\n",
    "    global pred_patch\n",
    "    \n",
    "    #If this is the first frame, initialize the view!\n",
    "    if i==0:\n",
    "        #Add all of our birds to the window\n",
    "        for flock in patches:\n",
    "            for patch in flock:\n",
    "                plt.gca().add_patch(patch)\n",
    "                npatches.append(patch)\n",
    "        #Add the predator to the window\n",
    "        plt.gca().add_patch(pred_patch)\n",
    "        npatches.append(pred_patch)\n",
    "        return npatches\n",
    "\n",
    "    #Calculate the wall avoidance force\n",
    "    wall_force = avoid_wall_force(flocks)\n",
    "    #Calculate the motive force that tries to keep the birds moving\n",
    "    mot_force = motive_force(flocks)\n",
    "    #Calculate the cohesion force\n",
    "    coh_force = cohesion_force(flocks)\n",
    "    #Calculate the alignment force\n",
    "    align_force = alignment_force(flocks)\n",
    "    #Calculate the separation force\n",
    "    sep_force = separation_force(flocks)\n",
    "    #Calculate the effect of the predator\n",
    "    pred_force, flee_force = predator_force(flocks, predator)\n",
    "    #Sum the forces to obtain the total force on each bird!\n",
    "    total_force = sum_arrays(wall_force,mot_force, coh_force, align_force, sep_force, flee_force)\n",
    "\n",
    "    #Now we update the position and velocity of the birds according to the force\n",
    "    updatePositionVelocity(flocks, total_force, patches, npatches)\n",
    "    #Do the same for the predator\n",
    "    pred_wall_force = avoid_wall_force([[predator]])[0][0]\n",
    "    updatePositionVelocityForPredator(predator, np.add(pred_wall_force,pred_force), pred_patch, npatches)\n",
    "\n",
    "\n",
    "    return npatches\n",
    "\n",
    "#Animate!\n",
    "anim = animation.FuncAnimation(fig, animate, \n",
    "                               init_func=init, \n",
    "                               frames=10000, \n",
    "                               interval=interval,\n",
    "                               blit=True)\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}