Boids: Revision Two

A friend of mine was in town this week, and haven't had much spare time to play with my side projects. When I finally got down to it I was able to redo my boids simulation to use operator overloading, boundary rules, and improved some of the rendering. I based my implementation off of a couple of different sources: a pseudo-code explanation of the original boids paper, a Python implementation, and this Java implementation. There is one remaining problem that I have yet to figure out. I originally used the rules laid out in the Python implementation, but the animation was very jerky when rule 3 took affect, and once the dots swarmed together they just sit their stationary. Finding no difference between his implementation and my own I started looking around for other examples that might give me an insight. That is when I found the previously mentioned Java implementation which resulted in a much more smooth trajectory. Unfortunately it still did not fix the problem where after swarming together the boids don't go anywhere. I'm going to continue to work on it but wanted to post my results thus far.

open System;;
open System.Windows.Forms;;
open Microsoft.FSharp.Collections;;

let BOID_COUNT = 40;;

let HEIGHT = 700.;;
let WIDTH = 1000.;;

let WALL = 50.;;
let WALL_FORCE = 30.;;

let MAX_VEL = 800.;;

let R1_CONST = 0.2;;
let R2_CONST = 0.2;;
let R2_RANGE = 100.;;
let R3_CONST = 0.2;;

let FPS = 24;;
let RECT_HW = 3;;


type vector = {
  x : double;
  y : double;
} with   

  member v.mag() = 
    (v.x ** 2.0) + (v.y ** 2.0) |> Math.Sqrt;

  static member ( + )(a:vector, b:vector) = 
    {x = a.x + b.x;
     y = a.y + b.y};

  static member ( - )(a:vector, b:vector) = 
    {x = a.x - b.x;
     y = a.y - b.y};

  []
  static member ( * ) (a:vector, b:int) = 
    let bd = Convert.ToDouble(b) in
    {x = a.x * bd;
     y = a.y * bd};

  []
  static member ( * ) (a:vector, b:double) = 
    {x = a.x * b;
     y = a.y * b};

  []
  static member ( ** ) (a:vector, b:int) = 
    let bd = Convert.ToDouble(b) in
    {x = a.x ** bd;
     y = a.y ** bd};

  []
  static member ( ** ) (a:vector, b:double) = 
    {x = a.x ** b;
     y = a.y ** b};

  []
  static member ( / ) (a:vector, b:int) = 
    let bd = Convert.ToDouble(b) in
    {x = a.x / bd;
     y = a.y / bd};

  []
  static member ( / ) (a:vector, b:double) = 
    {x = a.x / b;
     y = a.y / b};

  static member ( > ) (a:vector, b:double) = 
    a.mag() > b;

  static member ( < ) (a:vector, b:double) = 
    a.mag() < b;

end;;

type boid = {
  pos : vector;
  vel : vector;
} with   

  static member ( + )(a:boid, b:boid) = 
    {pos = a.pos + b.pos;
     vel = a.vel + b.vel};

  static member ( - )(a:boid, b:boid) = 
    {pos = a.pos - b.pos;
     vel = a.vel - b.vel};

  member b.distance b2 =
    Math.Sqrt(Math.Abs(((b.pos.x - b2.pos.x) ** 2.0) - ((b.pos.y - b2.pos.y) ** 2.0)));

  member b.update_pos() = 
    {b with
       pos = b.pos + b.vel};

  member b.neighbors(r, bl) = 
    [for n in bl -> 
       (b.distance(n), n)] |> 
        List.filter(fun (d:double, _) -> d < r) |> 
            List.sort(fun (d1:double, _) (d2:double, _) -> d1.CompareTo(d2)) |> 
                List.tl |> 
                    List.unzip;
    
end;;


let rand = new Random();;


let rand_double(max:double) = 
  Convert.ToDouble(rand.Next(Convert.ToInt32(max) - 1)) + rand.NextDouble();;


let vector_init = 
  {x  = 0.;
   y  = 0.};;


let boid_init = 
  {pos  = vector_init;
   vel  = vector_init};;


let boid_sum_list bl = 
  List.fold_left(+) boid_init bl;;


let boid_adjust_v_to_fps b = 
 { b with
      vel = b.vel / (Convert.ToDouble(FPS))};; 
 
  
let boid_limit_vel b = 
  let v = b.vel.mag() in
  match v with
    | _ when v > MAX_VEL ->
        let vl = v / MAX_VEL in
          { b with
              vel = b.vel / vl}

    | _ -> b;;



let boid_check_x_bounds b =
  match b with
    | _ when b.pos.x < WALL ->
        {b with 
           vel = {b.vel with 
                    x = b.vel.x + WALL_FORCE}};

    | _ when b.pos.x > (WIDTH - WALL) ->
        {b with 
           vel = {b.vel with 
                    x = b.vel.x - WALL_FORCE}};

    | _ -> b;;


let boid_check_y_bounds b =
  match b with
    | _ when b.pos.y < WALL ->
        {b with 
           vel = {b.vel with 
                    y = b.vel.y + WALL_FORCE}};

    | _ when b.pos.y > (WIDTH - WALL) ->
        {b with 
           vel = {b.vel with 
                    y = b.vel.y - WALL_FORCE}};

    | _ -> b;;


let boid_rule1 b_sum count =
  fun (b) ->
    {b with
       vel = b.vel + ((((b_sum.pos - b.pos) * (1. / (count - 1.))) - b.pos) * R1_CONST) };;


let boid_rule2 bl count = 
  fun (b:boid) ->
    let (_, nbl) = b.neighbors(R2_RANGE, bl) in
    let r2b = List.fold_left(fun b1 b2 -> b1 - (b2 - b)) boid_init nbl in
    { b with
        vel = b.vel + ((r2b.pos * (1. / (count - 1.))) * R2_CONST) };;


let boid_rule3 b_sum count =
  fun b ->
    {b with
       vel = b.vel + ((((b_sum.vel - b.vel) * (1. / (count - 1.))) - b.vel) * R3_CONST) };;


let boid_apply_rules b bl b_sum count = 
    b |> boid_check_y_bounds |> boid_check_x_bounds |> boid_rule1 b_sum count |> boid_rule2 bl count |>  
        boid_rule3 b_sum count |> boid_limit_vel |>  boid_adjust_v_to_fps |> fun b -> b.update_pos();;


let get_wall() = 
  let height = Convert.ToInt32(HEIGHT) in
  let width = Convert.ToInt32(WIDTH) in
  let wl = Convert.ToInt32(WALL) in
  [|new Drawing.Point(wl, wl);
    new Drawing.Point(wl, height - wl);
    new Drawing.Point(width - wl, height - wl);
    new Drawing.Point(width - wl, wl) 
  |];;


type BoidControl = class
    inherit UserControl as base
    
    val mutable boidList:(boid list);
    val brap:(Drawing.Point array);
    val count:double;
    val timer:Timer;
    
    new(bl:(boid list)) as this = {                                               
        boidList = bl;                                               
        count = bl |> List.length |> Convert.ToDouble;
        timer = new Timer();
        brap = get_wall();
    } then this.init();
    
    member this.init() = 
        let height = Convert.ToInt32(HEIGHT) in
        let width = Convert.ToInt32(WIDTH) in
        let size = new Drawing.Size(width, height) in

        base.Width <- width;
        base.Height <- height;
        base.MinimumSize <- size;
        base.MaximumSize <- size;
        base.BackColor <- Drawing.Color.White;
        base.SetStyle (ControlStyles.UserPaint, true);
        base.SetStyle (ControlStyles.DoubleBuffer, true);
        base.SetStyle (ControlStyles.AllPaintingInWmPaint, true);

        this.timer.Interval <- FPS;
        this.timer.Tick.Add(fun _ ->
                            this.tick(););

        this.timer.Start();

    member this.boids_update() =
        let b_sum = boid_sum_list this.boidList in
        let al = [for b in this.boidList -> async { return boid_apply_rules b this.boidList b_sum this.count }] |> Async.Parallel in 
        this.boidList <- (Async.Run(al) |> Array.to_list);
        ();
    
    override this.OnPaint e =
        let g = e.Graphics in
        List.iter(fun b -> 
            let x = Convert.ToInt32(b.pos.x) in
            let y = Convert.ToInt32(b.pos.y) in
            let bp1 = new Drawing.Point(x - RECT_HW, y - RECT_HW) in 
            let bp2 = new Drawing.Point(x - RECT_HW, y + RECT_HW) in 
            let bp3 = new Drawing.Point(x + RECT_HW, y + RECT_HW) in 
            let bp4 = new Drawing.Point(x + RECT_HW, y - RECT_HW) in 
            g.FillPolygon(Drawing.Brushes.Red, [|bp1; bp2; bp3; bp4|]);) this.boidList;
        
        g.DrawPolygon(Drawing.Pens.Blue, this.brap);
    
    member this.tick() =
        this.boids_update();
        this.Invalidate();
                
    
end;;


[]
do

  let height = Convert.ToInt32(HEIGHT) in
  let width = Convert.ToInt32(WIDTH) in
  let size = new Drawing.Size(width, height) in

  let bl = [for x in 1..BOID_COUNT -> {pos = {x = 50. + rand_double(WIDTH - 50.); 
                                              y = 50. + rand_double(HEIGHT - 50.)};
                                       vel = vector_init }] in
  let mw = new Form() in
  let bc = new BoidControl(bl) in
  mw.Controls.Add(bc);
  mw.ClientSize <- size;
  mw.MaximumSize <- mw.Size;
  mw.MinimumSize <- mw.Size;
  mw.Text <- "Boids";
  mw.Show();
  Application.Run(mw);;

Boids!

I've got revision one of my boids simulation done. I'll post more details and code later, but for now here is screen capture.

Asynchronous workflows in F#

Work has been really busy and I haven't had time to work on any of my side projects until last night, when I had chance to play with asynchronous workflows in F#. I've been wanting to explore to a greater extent what is out there in the world of concurrent and distributed computing, but in order to get started I needed to find some project that would be the catalyst. Not that having a project is absolutely necessary but it helps keep me on task and focused towards an end goal. Some of the ideas I've been tossing around were creating an engine to decrypt/encrypt PGP MIME messages, or build a better load/stress engine, both of which require a large up-front investment of time.

After about a week of going nowhere I stumbled upon the idea of doing swarm animations, and possibly throw in some additional AI behavior if I had time. The thing about doing swarm animations is that each entity in the swarm needs to move in relation to each of it's neighbors which involves a nearest neighbor search, one of those computationally costly and hard to tackle problems. I'm going to start with a really naive brute force approach and evolve it over time, into something more efficient. One of the things that I've been looking at is the Stolfi and Guibas method for creating Voronoi diagrams in order to speed up the nearest neighbor search. It's a divide an conquer algorithm so it should lend itself to parallelization, but that comes later, today I start with asynchronous workflows. The bit of code below had the following runtime on my Dell D620 dual-core laptop when executed using fsi:

Async time: 00:00:34.7309424
Sync time: 00:00:40.1761800

let largest = 800;;

type point = {
    x : double;
    y : double;
}

let rand = new Random();;

let rand_double(max) = 
  System.Convert.ToDouble(rand.Next(max - 1)) + rand.NextDouble();;

let rec quick_sort(x) =
  match x with
  | pivot :: rest ->
      let left, right = List.partition (( > ) pivot) rest in
      quick_sort left @ pivot :: quick_sort right
  | [] -> [];;

let distance(p1, p2) =
  Math.Sqrt(Math.Abs(((p1.x - p2.x) ** 2.0) - ((p1.y - p2.y) ** 2.0)));;

let time(f) =
  let start = DateTime.Now in
  f();
  DateTime.Now - start;;

let neighbors(p, pl) = 
  [for n in pl -> (distance(n, p), n)];;

let async_neighbors(p, pl) = 
  let a = Async.Parallel[for n in pl -> async {return (distance(n, p), n)}] in
  Async.Run(a);;

let sync_test() =
  let pl = [for x in 1..5000 -> {x=rand_double(largest);y=rand_double(largest)}] in
  [for p in pl -> neighbors(p, pl)] |> ignore;;

let async_test() = 
  let pl = [for x in 1..5000 -> {x=rand_double(largest);y=rand_double(largest)}] in
  let a = Async.Parallel[for p in pl -> async { return neighbors(p, pl) }] in
  Async.Run(a) |> ignore;;

let async_test2() = 
  let pl = [for x in 1..5000 -> {x=rand_double(largest);y=rand_double(largest)}] in
  let a = Async.Parallel[for p in pl -> async { return async_neighbors(p, pl) }] in
  Async.Run(a) |> ignore;;  

let cleanup() = 
  GC.Collect();
  GC.WaitForPendingFinalizers();;

do Console.WriteLine("Async time: " + any_to_string(time(async_test)));;  
do cleanup();;
do Console.WriteLine("Sync time: " + any_to_string(time(sync_test)));;

Calling asynchronous Twisted Python code from a blocking function

So I had to write a VAPI-XP test for Mercury/HP Quality Center, and the way that this particular test type works is that you write a script in Python, VBScript, or Perl. This script must have a TestMain function that is called and as soon as that function exits the test is over and everything exits. This test had to SSH into one of the machines in our virtualized test bed, and tell the machine to netboot in order to pick up the newest build of our server product.

I could have just used a call to os.popen() to some SSH client on our intermediary testing host but I had have bad experiences with that in the past, mostly weird behavior with OpenSSH and plink. I had written some SSH code with the Twisted module Conch before so I decided to just go ahead and use that. So got the SSH code up and running no problem, but how was I to block execution until the event was triggered letting me know that the session was complete? Furthermore how was I going to get my data back without doing something dirty like a global variable?

I talked to some guys on IRC and came up with the solution below. The SSH code itself isn't that interesting and I've only really included it for completeness the important bits to note is the call to "threads.blockingCallFromThread" in RunSSHCommand and note that I'm passing the deferred that I create in "_runSSHCommand" through each of the SSH classes (usually 'd' or 'self._d') and that it is finally being called in the method "CommandChannel.closed".

#!/usr/bin/python

class ClientTransport(transport.SSHClientTransport):

   def __init__(self, user, command, d):
       self._user = user
       self._command = command
       self._d = d

   def verifyHostKey(self, pubKey, fingerprint):
       #we don't care accept any host key
       return defer.succeed(1)

   def connectionSecure(self):
       #connection made, instiantiate the class that will handle authentication
       #and pass an instance of the class that embodies the user behavior as a param
       self.requestService(ClientUserAuth(self._user, ClientConnection(self._command, self._d)))


class ClientUserAuth(userauth.SSHUserAuthClient):

   def getPassword(self, prompt = None):
       # this says we won't do password authentication        
       return

   def getPublicKey(self):
       #return the public key which is defined up top as a string
       return keys.getPublicKeyString(data=SSH_PUB_KEY)

   def getPrivateKey(self):
       #return the pricate key which is also defined up top as a string
       return defer.succeed(keys.getPrivateKeyObject(data=SSH_PRIV_KEY))


class ClientConnection(connection.SSHConnection):

   def __init__(self, cmd, d, *args, **kwargs):
       connection.SSHConnection.__init__(self)
       self._command = cmd
       self._d = d


   def serviceStarted(self):
       self.openChannel(CommandChannel(self._command, self._d, conn=self))


class CommandChannel(channel.SSHChannel):
   name = 'session'

   def __init__(self, command, d, *args, **kwargs):
       channel.SSHChannel.__init__(self, *args, **kwargs)
       self.command = command
       self.d = d
       self.data = ""

   def channelOpen(self, data):
       #send an execute request passing the user supplied string.
       self.conn.sendRequest(self,
                             'exec',
                             common.NS(self.command),
                             wantReply=True
                             ).addCallback(self._gotResponse)

   def _gotResponse(self, _):
       #command has returned, send an EOF in order to terminate the connection
       self.conn.sendEOF(self)

   def dataReceived(self, data):
       #append the output data to the string
       self.data = self.data + data

   def closed(self):
       #connection closed, execute the callback on the deferred
       #that we have been passing around
       self.d.callback(self.data)


class ClientCommandFactory(protocol.ClientFactory):
   #factory class for our SSH protocol
  
   def __init__(self, user, command, d):
       self._user = user
       self._command = command
       self._d = d

   def buildProtocol(self, addr):
       protocol = ClientTransport(self._user, self._command, self._d)
       return protocol


def RunSSHCommand(hostname, port, user, cmd):
   #Run an SSH command as a user on a given host, we only use key authentication
   thread.start_new_thread(reactor.run, (False,))
   try:
       return threads.blockingCallFromThread(reactor, _runSSHCommand, hostname, port, user, cmd)
   except:
       return None
  
   reactor.stop()

def _runSSHCommand(hostname, port, user, cmd):
   d = defer.Deferred()
   factory = ClientCommandFactory(user, cmd, d)   
   reactor.connectTCP(hostname, 22, factory)
   return d
  

if __name__ == "__main__":
   RunSSHCommand(hostname, port, user, "register_netboot")