ofile := "/tmp/toto" ;;

record pre_boid = { x:float, y:float, t:float, w:int } ;;
record boid = pre_boid + { w = 1 } ;;
geoprox population((20, 20), 5.0) = fun b -> (b.x, b.y) ;;

fun cons_random_boid(id, pop) = { x = random(100.), y = random(100.), t = random(2.0 * PI), w = 1 } :: pop ;;

S0 := fold(cons_random_boid, population:(), 100) ;;

fun dist2(b1,b2) = (
  let dx = b1.x - b2.x
  and dy = b1.y - b2.y
  in
    dx * dx + dy * dy
) ;;

fun add_boids(b1,b2) =
  { x = b1.x + b2.x,
    y = b1.y + b2.y,
    t = b1.t + b2.t,
    w = b1.w + b2.w }
;;

fun move_boid[dt=0.1](b) =
  let x = (b.x + dt * cos(b.t) + 100.0) % 100.0
  and y = (b.y + dt * sin(b.t) + 100.0) % 100.0 in
    b + { x = x, y = y }
;;

fun too_close(b1, b2) = dist2(b1,b2) <= 1.0 ;;

fun too_far(b1, b2) = dist2(b1,b2) >= 16.0 ;;

trans rules = {

  b / neighbors_exists(too_close(b), b) => (
    let g = neighbors_fold(add_boids, { x = 0., y = 0., t = 0., w = 0 }, b) in
    let dx = b.x - g.x / g.w
    and dy = b.y - g.y / g.w  in
    let t = if dy == 0. && dx <= 0. then PI else 2. * atan(dy / (dx + sqrt(dx*dx + dy*dy))) fi in
    let b' = b + { t = t } in
      !! g.w > 0;
      move_boid(b')
  );

  b / neighbors_forall(too_far(b), b) => (
    let g = neighbors_fold(add_boids, { x = 0., y = 0., t = 0., w = 0 }, b) in
      if g.w == 0 then move_boid(b)
      else
        let dx = g.x / g.w - b.x
        and dy = g.y / g.w - b.y in
        let t = if dy == 0. && dx <= 0. then PI else 2. * atan(dy / (dx + sqrt(dx*dx + dy*dy))) fi in
        let b' = b + { t = t } in
          move_boid(b')
      fi
  );

  b => (
    let g = neighbors_fold(add_boids, b, b) in
    let b' = b + { t = g.t / g.w } in
      move_boid(b')
  );

} ;;

fun output(S) = (
  if (iteration % 10 == 0)
  then
    ofile << " { " ;
    foreach b in S do (
      ofile << " { " << b.x << " " << b.y << " 0 0 0 } ";
      ofile << " [ " << b.x << " " << b.y << " " << (b.x + 5*cos(b.t)) << " " <<(b.y + 5*sin(b.t)) << " 1 0 0 ] "
    );
    ofile << " }\n"
  fi;
  S
) ;;

rules[iter = 10000,
      prelude = output,
      interlude = output,
      postlude = output
     ](S0) ;;