Sunday, August 28, 2011

firing and forgetting... Akka continuations...

Strange title. I found myself with my hands full for a lot of time,and dumbly had nothing to talk about. Reading Scala In Action, Clojure In Action, and reading the Haskell books on the table will take a lot of time.Don't forget I am getting crazier about Lisp. So what ? I found myself worried because I dropped what I had started reading about Akka.
Like many people, if I do not practice regularly I have to rollback again to the very beginning of the study and try again. Unfortunate. This time I chose to come back to the basic and bought Gul Agha book about the Actor Model paradigm. This will help me in parallel to acquire some inherent knowledge about the paradigm while reaching the last chapter of Scala in action, about fault tolerant systems and all the rest. This is my purpose with Akka actors. Today's purpose is a "back-to-the-basics" reminder.

If by any chance Jonas Boner reads this article , may he forgive me for publishing a so basic level example, Akka deserves the best. Starting again on Akka, I downloaded the 1.2RC3 version and decided in order to learn my way again to implement one of the examples of the chapter 3. The famous recursive factorial example. Although this is not exactly the subject of my article, as a reminder, this is how a recursive factorial can be processed from an actor point of view:



The actor 1 receives a message containing the value of the number n to be processed.
On each message reception:
  • it does create (dashed arrow) some new actor that will accumulate a number k
  • decrements the received number and resend to itself a message if the new number is not 1
  • otherwise sends a process message to the last created accumulated actor which is storing 1

Accumulating actors will expect for a processing message containing a number to be multiplied by the number they are actually storing. They will propagate the result to the accumulating actor containing
the next number or send back the result. So In the case you would want to process factorialOf(3), then the accumulating actor 2 would store 3, the actor 3 would store 2 and the fourth would store 1.

When actor 1 fires a message to actor 4 in order to complete the operation then the whole chain is processed and the result is sent to actor 6.

I willingly implemented a fork solution splitting the task.

I loved the idea although the performances are very poor. I loved it because this model of chaining actors reflects the genuine Alan Kay idea of object oriented message passing. Thanks to Akka you really create instances of objects, each owning a single responsibility. That is Object Oriented programming a la Smalltalk and not class oriented programming a la "whatever you want" (please don't misunderstand me, I love Java too).

I did not like the performances to be poor (5 times a tail recursive version) so I decided to reduce the overload of message passing reducing the number of accumulating actors. The bigger your number to process, the bigger your number of accumulating actors, so your process time grows exponentially with the input number value. I created a fork of the tasks to process:


The input actor 1 branches the looping actors 2, 4, 6. The looping actors creates their storing actors and at the end of the process all the work is sent to the continuing actor 8. Its purpose being to synchronize all the jobs, as being also a mean of control afterwards, this actor for me is nothing more than a continuation. This observation, lead me to consider a 1977 article of Carl Hewitt, where this taxonomy is used. Still reading this fascinating article.

In order to challenge the system I created a first set of tests, and for comparing purpose, implemented both a tail recursion solution and the actor solution. Here come the tests:

import org.specs2._
import MathBox._
import akka.actor.Actors
import org.perf4j.log4j.Log4JStopWatch

final class FactorialSpecification extends Specification {
  def is =
    "Factorial Specification" ^
      p ^
      "Tail recursive factorial should" ^
      "provide the expected result !!! " ! e1 ^
      p ^
      "Acting Factorial should" ^
      "be processed equal to tail recursive one " ! AssertOnActors(2048).e1


  def e1 = factorialOf(5).should(beEqualTo(120))

  case class AssertOnActors(number: Int) extends specification.After {
    def e1 = {
      val stopWatch = new Log4JStopWatch();
      Range(0, 1000).foreach {
        value =>
          actingFactorialOf(number)
          stopWatch.lap("parallel");
          factorialOf(number)
          stopWatch.lap("serial");
      }
      actingFactorialOf(number).should(beEqualTo(factorialOf(number)))
    }

    def after {
      Actors.registry().shutdownAll()
    }
  }

}


And yes although I do not like frameworks, because it makes you dumb, I used perf4j to measure performances (so I am a little dumb). The assertions validate the processing of the tail recursive factorial and, then validate the actors solution thanks to the tail recursive one.

The entry points of the processing are located into the MathBox object:

import annotation.tailrec

object MathBox {
  def actingFactorialOf(number: Int): BigInt = {
    val factorial = FactorialActor()
    (factorial ? FactorialOf(number)).get.asInstanceOf[FactorialResult] match {
      case FactorialResult(result) =>   result
      case _ => -1
    }
  }


  def factorialOf(number: BigInt): BigInt = {
    @tailrec
    def recFactorialOf(accumulated: BigInt, counter: BigInt): BigInt =  {
      if (counter <= 1) accumulated
      else recFactorialOf(accumulated * counter, counter - 1)
    }

    recFactorialOf(1, number)
  }
}

Using the Scala embedded @tailrec annotation verifies that the method will be compiled with tail call optimization, and if not will throw some exception. As I did not want to clutter the tests with timing I used the new symbol provided in the 1.2-RC3 to get a future after sending a message: the '?'. What I do send is a FactorialOf(number) message while expecting a FactorialResult(result). Nothing terrific.
I told you this article was easy, but I had to write something blah blah blah... :)
 So comes the solution:

import akka.actor.Actor._
import akka.actor.{PoisonPill, ActorRef, UntypedChannel, Actor}
import scala.Option

case class FactorialOf(number: Int)

case class FactorialResult(result: BigInt)

case class StoreResult(number: Int, channel: Option[UntypedChannel])

case class Cumulate(number: Int)

case class WorkOn(min: Int, max:Int, next: UntypedChannel)


class AccumulatingActor(val target: UntypedChannel) extends Actor {
  var stored: BigInt = 1

  protected def receive = {
    case Cumulate (operand)=> stored = stored * operand
    case FactorialResult(operand) =>
      target ! FactorialResult(stored * operand)
      self.getChannel ! PoisonPill
    case _ =>
  }
}
object AccumulatingActor {
  def apply(target: UntypedChannel): ActorRef = {
    actorOf(new AccumulatingActor(target)).start()
  }
}

class LoopingActor extends Actor {
  var min: Int = _
  var intermediate: UntypedChannel = _

  protected def receive = {
    case WorkOn(lower, max, target) =>
      min = lower
      intermediate = AccumulatingActor(target)
      become(decrementer)
      self ! FactorialOf(max)
  }

  protected def decrementer: Receive = {
    case FactorialOf(number) =>  {
      intermediate ! Cumulate(number)
      if (number > min) {
        self ! FactorialOf(number - 1)
      } else {
        intermediate ! FactorialResult(1)
     }
    }
  }
}

object LoopingActor {
  def apply() = actorOf[LoopingActor].start()
}

class ContinuingActor(val expectedMessages: Int,val  target: UntypedChannel) extends Actor {
  var counter = expectedMessages
  var result: BigInt = 1

  protected def receive = {
    case FactorialResult(value) =>
      counter = counter - 1
      result = result * value
      if (counter == 0) {
        target ! FactorialResult(result)
        self ! PoisonPill
      }
  }
}

object ContinuingActor {
  def apply(expectedMessages: Int,target: UntypedChannel) = {
    actorOf( new ContinuingActor(expectedMessages, target)).start()
  }
}

class FactorialActor extends Actor {
  def coreNumber = Runtime.getRuntime.availableProcessors() * 8;

  def dispatch(number: Int, bunch: Int, channel: UntypedChannel) {
    bunch match  {
      case 0 =>
        LoopingActor() ! WorkOn(1, number, Some(ContinuingActor(1, channel)))
      case _ =>
        val continuation = ContinuingActor(coreNumber, channel)
        val lastBoundary = Range(0, coreNumber - 1).foldLeft(0){ (acc, step) =>
          LoopingActor() ! WorkOn(acc + 1, acc + bunch , continuation)
          acc + bunch
        } + 1

        LoopingActor() ! WorkOn(lastBoundary, number, continuation)
    }
  }

  def asBunch(forMax: Int): Int = forMax / coreNumber

  protected def receive = {
    case FactorialOf(number) => dispatch(number, asBunch(number), self.getChannel)
  }
}

object FactorialActor {
  def apply() = {
    actorOf[FactorialActor].start()
  }
}

Waoooooooh !! Too long. Don't panic. There is no challenge in it. Let us start with the FactorialActor actor !

class FactorialActor extends Actor {
  def coreNumber = Runtime.getRuntime.availableProcessors() * 8;

  def dispatch(number: Int, bunch: Int, channel: UntypedChannel) {
    bunch match  {
      case 0 =>
        LoopingActor() ! WorkOn(1, number, Some(ContinuingActor(1, channel)))
      case _ =>
        val continuation = ContinuingActor(coreNumber, channel)
        val lastBoundary = Range(0, coreNumber - 1).foldLeft(0){ (acc, step) =>
          LoopingActor() ! WorkOn(acc + 1, acc + bunch , continuation)
          acc + bunch
        } + 1

        LoopingActor() ! WorkOn(lastBoundary, number, continuation)
    }
  }

  def asBunch(forMax: Int): Int = forMax / coreNumber

  protected def receive = {
    case FactorialOf(number) => dispatch(number, asBunch(number), self.getChannel)
  }
}

object FactorialActor {
  def apply() = {
    actorOf[FactorialActor].start()
  }
}

There, on reception of a FactorialOf(number) message:

protected def receive = {
    case FactorialOf(number) => dispatch(number, asBunch(number), self.getChannel)
  }

The FactorialActor
  • creates the continuation
  • fork the looping actors 
The FactorialActor grossly cut the number into bunches or range of number to be processed using the raw heuristic:

def coreNumber = Runtime.getRuntime.availableProcessors() * 8;

As all of the actors in this example, the FactorialActor instance is created using a companion object. The created looping actors are started using a WorkOn message:

LoopingActor() ! WorkOn(start, end , continuation)

holding the inclusive boundaries of the range to work on, and the continuation to be passed the results to. So what about the continuation ?

class ContinuingActor(val expectedMessages: Int,val  target: UntypedChannel) extends Actor {
  var counter = expectedMessages
  var result: BigInt = 1

  protected def receive = {
    case FactorialResult(value) =>
      counter = counter - 1
      result = result * value
      if (counter == 0) {
        target ! FactorialResult(result)
        self ! PoisonPill
      }
  }
}

object ContinuingActor {
  def apply(expectedMessages: Int,target: UntypedChannel) = {
    actorOf( new ContinuingActor(expectedMessages, target)).start()
  }
}

As a good kamikaze actor, the continuing actor will store the incoming partial results, finalizing the multiplication in order to send it to its target, before committing suicide with a poison pill.
Therefore the continuing actor also plays the part of a forwarding actor sending back the result to the client, hiding its reference form the other actors view.
We never have to send back the result to previous actor like we do in standard fork/join frameworks because we have the control as a continuation.
So, we avoid cluttering the code of the pattern matching in the emitting previous actors preventing them from handling a return result.

 Remain two actors more. First the looping actor in charge of looping into the range boundaries:

class LoopingActor extends Actor {
  var min: Int = _
  var intermediate: UntypedChannel = _

  protected def receive = {
    case WorkOn(lower, max, target) =>
      min = lower
      intermediate = AccumulatingActor(target)
      become(decrementer)
      self ! FactorialOf(max)
  }

  protected def decrementer: Receive = {
    case FactorialOf(number) =>  {
      intermediate ! Cumulate(number)
      if (number > min) {
        self ! FactorialOf(number - 1)
      } else {
        intermediate ! FactorialResult(1)
     }
    }
  }
}

object LoopingActor {
  def apply() = actorOf[LoopingActor].start()
}

I chose to apply a behavior hot swap, splitting the behavior of being initialized from the behavior of self looping on decremented FactorialOf(number). The hot swap can be managed using the nice become method switching to the decrementer behavior.
During initialization, a freshly created AccumulatingActor instance is passed the continuation and the swap is achieved The decrementer behavior self loop until the lower boundary is reached , then the looping actor fires a message of end of calculus using a FactorialResult(1) message.
A FactorialResult(1) message flags the end of a processing for a range. One should notice that in his 1986 book, Gul Agha already exposes the change of behavior presenting some prototype language which syntax is very close to Erlang or Akka syntaxes. On the Accumulating actor side:

class AccumulatingActor(val target: UntypedChannel) extends Actor {
  var stored: BigInt = 1

  protected def receive = {
    case Cumulate (operand)=> stored = stored * operand
    case FactorialResult(operand) =>
      target ! FactorialResult(stored * operand)
      self.getChannel ! PoisonPill
    case _ =>
  }
}
object AccumulatingActor {
  def apply(target: UntypedChannel): ActorRef = {
    actorOf(new AccumulatingActor(target)).start()
  }
}

the accumulating actor simply processes a multiplication and forward the result to the continuation actor. 
We are done. Execution, tests green (joking, it took me several shots ;)) A sample of average statistics

Performance Statistics
Tag          Avg(ms)         Min         Max     Std Dev       Count
parallel     3,6             3            19         1,0         485
serial       12,1            11           73         2,9         485

Performance Statistics  
Tag          Avg(ms)         Min         Max     Std Dev       Count
parallel     10,1           3           1082        77,8         191
serial       12,2           12            18         0,7         191

Performance 
Tag          Avg(ms)         Min         Max     Std Dev       Count
parallel     3,7              3           6          0,8         588
serial       12,2            11          66          2,3         587


exhibits warm up phases, than in nominal conditions a nice ratio from 1/4 concerning the increase of performances. Note that the tail recursion performances (serial tag) are very close to the loop performances I benched using Java.

I realised experimentations with the ForkJoin Framework, as I was using JDK7 in order to run Scala (2.9.1RC4). Performances are really close in favor to the ForkJoin framework, but one should notice that I did not try to bench greater values and did not respect a very serious protocol.

The nice thing, as far as I am concerned, is the ability to adopt a perfectly understandable object oriented approach with the actor paradigm. We used CPS as in functional programming languages, we have hidden the client channel, we have separated concerns etc...

Don't know what the end of the book will reveal but can't wait to read it. Hoping the next time we will talk about actors will be about fault tolerance, have a nice day.

Be seeing you :)!!!!




Monday, August 15, 2011

Sipping some Monte-Carlo in Scala and Clojure

Howdy buddies, a new one from the SICP. The lecture's nice and suits well the learning curve of all kinds of functional programming languages. Naturally that does concern Clojure, as a Lisp-1 natural language, and Scala. I recently watched the lecture concerning the (very) few benefits of assignments. There, Gerald Jay Sussman presented a simple example of a Monte Carlo simulation dedicated to the processing of the number Pi.

First, thanks a lot to Debasish Gosh who kindly accepted to answer my questions and inspect my code too. May he be assured of all my gratitude. Exchanging ideas with a so skilled person like him was quite an experience and enlightened these times of trouble for me. I hope I have reported all the modifications he suggested, may he forgive me if I did not. I may not have found yet some Master Craftsman in Europe in order to teach me about actors and Functional programming (40 year old may be too...old), but the kind help of this person was a real support.

Back to Pi. Some approximated value of pi can be quickly computed if based on the Cesaro demonstration that the square of Pi is inversely proportional to probability that two integers chosen at random will have no factors in common. More clearly:

P(gcd(N1,N2) == 1) = 6 /  π *  π

Because if so, their gcd is 1. Quoting the Wikipedia article "Monte Carlo methods (or Monte Carlo experiments) are a class of computational algorithms that rely on repeated random sampling to compute their results". So running a Monte Carlo experiment in order to process the Pi number can be resumed in running a gcd calculus on random natural integer values till our result converges to some expected value.
Therefore, random number generation takes place here, if we consider running a gcd function on a series of nearly random integer. In order not to make things too complex, let's assume that the action of generating numbers can be described as the process of creating a string of numbers extrapolated from a suite taking its roots from a seed value. I picture that as a mathematical suite :

Vn = f(Vn-1) so

Vn = f(f(f(.................f(Vo)

where Vo = some_seed

where at each step one need to cache the new generated value in order for the following one to be processed and so on. (Ok my mathematical symbolic vocabulary is limited but I do not want to copy n'paste Wikipedia)
This is a typical problem solved with the help of assignable variables. Of course some solution does exist which does not depends on some external variable but the result is much more cluttered (have a look there).

Having installed Leinengen (version 1.6.1 is very nice with repl at its top), I now have all the tools I need to challenge my code with tests. In the Clojure package cesaro.test I created a suite of small tests in a core_spec.clj file. Starting content is :

(ns cesaro.test.core-spec
  (:use cesaro.core)
  (:use clojure.test))

Basic. As in the cesaro.test package, in the core_spec.clj file content, my namespace is cesaro.test.core-spec. I claim there my intent to use the tools of the clojure.test namespace in order to challenge the functions inplemented in the cesaro.core namespace (probably the content of a core.clj file in a cesaro package... got it? :))

What I need is a working gcd function first. So are the very dumb tests used to create it:

(deftest gcd-with-matching-known-numbers-should-return-value
	(is (= (gcd 1989 867) 51)))

(deftest gcd-with-matching-some-other-numbers-should-return-value
	(is (= (gcd 36 27) 9)))

(deftest gcd-with-prime-numbers-should-return-one
	(is (= (gcd 23 11) 1)))

(deftest gcd-with-unordered-numbers-should-return-gcd
	(is (= (gcd 11 23) 1)))

The last test appeared later as I used some Euclide method so to process the gcd and had to face slight problems due to my expectation of ordered parameters :) (told you I was dumb...but still learning).
This leads me to:

(defn euclide-gcd [numberOne numberTwo]
	(cond 
		(= numberTwo numberOne) numberOne
		(= numberTwo 1) 1
		:else (let [dif (- numberOne numberTwo)]
			(if (< dif numberTwo)
				(recur numberTwo dif)
				(recur dif numberTwo)))))

(defn gcd [x y]
	(if (< x y) (euclide-gcd y x) (euclide-gcd x y)))

The euclide-gcd function expects ordered parameters. I do not like the conception and what I hate most are the cluttered cond/if branches. This is something I intend to avoid in the future. I hate branches as they do remind me of goto. Clojure deserves better than that.
Then I need some random generator. A mechanism that would allow me to generate a seed and the derived random numbers. The tests look like:

(deftest seeds-with-two-invocations-should-differ
	(is (not= (seed) (seed))))

(deftest random-with-two-invocations-should-differ
	(is (not= (random) (random))))

The implementation I provided :

(defn seed [] (rand-int Integer/MAX_VALUE))

(defn marsaglia-first[value]
	(+ (* 36969 (bit-and value 65535)) (bit-shift-right value 16)))

(defn marsaglia-second[value]
	(+ (* 18000 (bit-and value 65535)) (bit-shift-right value 16)))

(defn marsaglia-sum [x y] (+ (bit-shift-left x 16) y))

(def random 
	(let [x (ref (seed)) y (ref (seed))]
		(defn update []
			(dosync
				(alter x marsaglia-first)
				(alter y marsaglia-second)
				(marsaglia-sum @x @y)))))

And yes, I adopted a Marsaglia algorithm (see wikipedia previous reference). Although quite idiomatic, I find the result not as elegant as the Scheme solution where the set! form allows for the modification of a declared variable in a let form expression:

(define rand
  (let ((x random-init))
    (lambda ()
      (set! x (rand-update x))
      x)))

The Scheme code shape appears to be more concise. This concision is not due to the fact that I chose a Marsaglia algorithm involving two parameters instead of one, but it finds its origin in the fact that - by nature -, variables set by a let form in Clojure are immutable. The only solution I found was to use the idiomatic form of the references in Clojure, in order to both alter the references wrapped values. Any suggestion will be welcomed.

Frustrated by the two previous pieces of code, I decided to force myself to write a Monte carlo small running function, without no branching at all. A dumb test is to provide my upcoming function with always/never failing tests in order to check its limits:

(deftest montecarlo-with-always-successfull-simulation-should-return-1
	(defn always-successfull [] true)
	(is (= 1 (montecarlo 10 always-successfull))))

(deftest montecarlo-with-never-successfull-simulation-should-return-0
	(defn never-successfull [] false)
	(is (= 0 (montecarlo 10 never-successfull))))


The montecarlo function accepts two input parameters, the first being the number of attempts and the second the simulation to be applied. The returned result is the ratio of passed simulation versus the whole number of simulations.

No branches so? I would lie saying that that one was a piece of cake for a beginner in functional programming as I am. I had to switch into Scala, than come back to Clojure with this following precise implementation:

(defn with-result-updates []
	(let [increments {true [1 0] false [0 1]}]
		(fn [status] (increments status))))

(defn run-montecarlo [updates in-simulation]
	(fn [[passed missed] number-of-essays]
		(vec (map + (updates (in-simulation)) [passed missed]))))

(defn montecarlo [essays in-simulation]
	(let [runs (run-montecarlo (with-result-updates) in-simulation) from-start [0 0]]
		(/ (first (reduce runs from-start (range essays))) essays)))

Of course the entry point is the montecarlo function. As I know the number of essays to run, all I need is to iterate and not recur over a range of essays:

(range essays)

Meanwhile, at each step, I can run a simulation (the purpose of the run-montecarlo function) that will update a vector of statistics, provided as [0 0] at the very beginning:

from-start [0 0]

The first element being the number of passed essays and the second the number of failed ones. The reduce form (the equivalent of the Scala foldLeft) aims to produce a vector of statistics. With reduce, one can produce whatever he wants, even another list, from the driving list. The purpose of the with-result-updates function is to create a lambda function instance (so a procedure) capable of producing a vector of the deltas to be applied whether a simulation has succeeded or not:
  • [1 0] on success
  • [0 1] on failure
One instance of the procedure is created once, so is the hashmap embedding the possible results:

{true [1 0] false [0 1]}

We close onto one immutable single instance of the hasmap, embedded into the frame (context of execution) of the created procedure. The principle of closing applies too in the montecarlo function where we both close onto the previous described procedure instance and the simulation to be applied:

(let [runs (run-montecarlo (with-result-updates) in-simulation) from-start [0 0]]

The final test to run to check on pi can be :

(deftest square-pi-with-enough-essays-must-be-close-to-9-9
	(let [result (square-pi)]
		(println "square-pi" result)
		(is (< (- result 9.9) 0.025))))

where we assert on the vicinity of the values of square pi and 9.9. Come then a natural implementation:

(defn cesaro-test []
	(let [value  (gcd (random) (random))]
		(= 1 value)))

(defn square-pi []
	(/ 6 (montecarlo 16192 cesaro-test)))

Works nice. What about Scala ? Scala helped me to find the no-branch version of the montecarlo method. For the trained Java eye, Scala is a pleasant bridge to take in order to embrace good functional programming habits to be adopted in any other functional programming language. Don't misunderstand me. Clojure and Scheme overwhelm me with thrilling sensations each time I do practice them. They also show me I was rambling in the dark before.
So, going back again to Scala (did I say I bought the Scala T-shirts and Teddy bear? ), I first had to write Random number generator tests. This was an opportunity to try Specs2:

import org.specs2._
import RNG._

final class RNGTest extends Specification{ def is =

  "RNG specification"               ^
                                    p^
  "Seed Generator should"           ^
  "Generate two different seeds"    ! e1^
  "Generate positive numbers"       ! e2^
                                    p^
  "Number Generator should"         ^
  "Generate two different numbers"  ! e3


  def e1 = seed should (not (beEqualTo(seed)))

  def e2 = {
    val rand = RNG()
    Range(0, 100).map((value: Int) => rand() should (beGreaterThan(0)))
  }

  def e3 = {
    val rand = RNG()
    rand() should not (beEqualTo (rand()))
  }

}

leading me to :

import util.Random
import  Long._
import scala.math._
class RNG(var seedOne: Int, var seedTwo: Int ) {

  def apply(): Int = {
    seedOne = 36969 * (seedOne  & 65535) + (seedOne  >> 16)
    seedTwo = 18000 * (seedTwo & 65535) + (seedTwo >> 16);
    abs((seedOne << 16) + seedTwo)
  }
}

object RNG {
  Random.setSeed(MaxValue)


  def seed : Int = {
    abs(Random.nextInt())
  }

  def apply(): RNG = {
    new RNG(seed, seed)
  }
}

where the seed is provided by the native number generator. I also need a GCD so let's test it

import org.specs2._
import com.promindis.montecarlo.MathModule._

final class GCDTest extends Specification {
  def is =
    "GCD calculus specification"                          ^
                                                          p^
      "GCD with macthing numbers should"                  ^
      "Find a first known result"                         ! e1 ^
      "Find a second known result"                        ! e2 ^
      "Find known result with reversed paramteres"        ! e3 ^
                                                          p^
      "GCD with non macthing numbers should"              ^
      "assert on known mismacth"                          ! e4 ^
      "assert on known mismacth with invert parameters"   ! e5


  def e1 = gcdOf(1989, 867) should be equalTo (51)

  def e2 = gcdOf(36, 27) should be equalTo (9)

  def e3 = gcdOf(27, 36) should be equalTo (9)

  def e4 = gcdOf(23, 11) should be equalTo (1)

  def e5 = gcdOf(11, 23) should be equalTo (1)
}

and then the solution:

object MathModule {

 //.........

  def sorted(firstNumber: Int, secondNumber: Int): (Int, Int) = {
    if (firstNumber < secondNumber) (secondNumber, firstNumber) else (firstNumber, secondNumber)
  }

  def gcdOf(firstNumber: Int, secondNumber: Int): Int = {
    def gcdOf(pair: (Int, Int)): Int = {
      pair match {
        case (y, x) if (x == y) => x
        case (y, x) if x == 1 => 1
        case (y, x) if x < y => gcdOf(sorted(y - x, x))
        case (y, x) if x < y => gcdOf(sorted(x - y, y))
      }
    }

    gcdOf(sorted(firstNumber, secondNumber))
  }
}


One can admire the beauty of the self expressive pattern matching in Scala. Finally I need a Monte Carlo simulator. Writing the tests again with the beautiful specs2:

import org.specs2._

final class MontecarloTest extends  Specification{ def is =

  "Montecarlo simulation dshould"                     ^
                                                      p^
  "have 100% success with always successful scenario" !e1 ^
  "have 0% success with always failing scenario"      !e2

  def e1 = {
    def alwaysSuccessful(): Boolean = true
    val stats = Montecarlo.simulation(alwaysSuccessful, 100)
    stats._1 should (beEqualTo(100))
    stats._2 should (beEqualTo(stats._1))
  }

  def e2 = {
    def neverSuccessful(): Boolean = false
    val stats = Montecarlo.simulation(neverSuccessful, 100)

    stats._3 should (beEqualTo(stats._1))
  }

}

helped me to get to :

object Montecarlo {
  type test = () => Boolean

  val update: Map[Boolean, List[Int]] =
    Map[Boolean, List[Int]](true -> List(1,0), false -> (List(0, 1)))

  def simulation(onRunning: test, essays: Int): (Int, Int, Int) = {
    val results: List[Int] = Range(0, essays).foldLeft(List(0, 0)) {
      (stats: List[Int], index: Int) =>
        (stats, update(onRunning())).zipped.map[Int, List[Int]](_ + _)
    }
    (essays, results(0), results(1))
  }
}

Here the result is provided as 3-Tuple, returning the number of essays, the number of passed essays , then the number of failed essays. I suffered only on the map method application after zipping the resulting lists. The compiler seemed in need of some help with explicit typing in order to infer the returned type. The astute reader will note that we used the same trick as in Clojure, storing the increments values definitions for the two success/failure scenarii into a Map. I have all the tools I need to run a Cesaro test:

import org.specs2.mutable.Specification

final class PiResolutionTest extends Specification{

  "Resolution of Pi" should {
    "be close to 9.9 " in {
      val epsilon = (9.9 - MathModule.estimatPiSquare())
       epsilon should(beLessThan(0.05))
    }
  }

}

and challenge it:

object MathModule {

  def forCesaro(random: RNG): () => Boolean = {
    () =>  gcdOf(random(), random())  == 1
  }


  def estimatPiSquare(): Double = {
    val stats = Montecarlo.simulation(forCesaro(RNG()), 10000)
    println(stats)
    val ratio: Double = int2double(stats._2) / stats._1.asInstanceOf[Double]
    println(ratio)
    println(6.0 / ratio)
    (6.0 / ratio)
  }
//..............

}

What was learnt ?

Well, living without branches is not easy, but worthwhile because it helps in learning how to use and reuse functional programming bricks. But I also learnt I have to digg into the RNG code so to find a smarter way to generate my numbers in Clojure. Maybe a stream oriented approach would be better...
Nice, got to finish chapter 11 of the Joy of Clojure and read one chapter more in Gul Agha's Actors book.


Be seeing you ! :)


Saturday, August 6, 2011

Scala tribute to SICP, Paul Henderson and... George.

Houdy, strange title. Took me time to produce something acceptable. Although not satisified of the code I decided to write something about some illustrative micro DSL I have been produced while watching some SICP lecture .As usual let us start with something completly different. Thanks to Mario Fusco and Jonas Boner for referencing my Hanoï tower article. This provided me with my fifteen minutes of fame, attracting a number of visitors very close to 5k last month.

I was quite shocked when I learnt that SICP was under attack,as I started watching the online videos a few months ago. I had the feeling of discovering again my job. Progressing step by step in the exploration of functional programming, I try to reproduce each difficult example when possible in Scala and Clojure. The Hanoi tower article came from there. Now I will take again the book and read it from start to end. The best way of protesting. I know I promised Debasish Ghosh to read the Haskell manual he suggested. I will,  it will only double my load of work, that's all. But SICP matters !!

The best medecine to cure this shock state was to provide another example that I wanted to implement for long. In lecture 3A, Hal Abelson presents some implementation of Peter Henderson's article presenting a functional geometry language allowing him to denote the structure of an Escher woodcut. Nice article, I did not understand everything. But at least I am progressing. I may not go as far as Hal Abelson, but I wanted to make a try in Scala.

The idea is to conceive some easy algebric language allowing us to express the shape of a graphic representation rather than detailing how it is built. This expression adopts the shape of an embedded language into Scala (or Clojure/Lisp etc...)
Functional programming comes to help. This is easily done if your elementary bricks are not modeled as objects but as procedures (runtime constructions of lambdas).

Imagine one moment that the most basic element of your construction would be a picture. Not an object per se as perceived in the OOP world, but the ability to proceed a rectangle zone in some graphic area (Java2D developers may see where I am going to).
The lambda definition hides the craft of drawing a category of drawings. At some point in time, a procedure will be created from the lambda, closing onto some basic draw description. Now your application knows how to draw a specific picture, whatever is your target rectangle area. This knowledge is applicable as many times as you want.

Creating well named combining functions (combinators)- in essence higher order functions manipulating other functions- one can setup an embedded language allowing for a clear expression of its intents as a combination  of higher order functions.

This is where George comes to play. Well, Georges is the graphical example provided by Harold Abelson:


In the algebraic language the following expression

beside(picture, above(empty(), picture, byHalf), byHalf)

would draw some provided pictures,  side by side in an input rectangle area, allocating half of the rectangle width to each of the picture version. The right picture is reduce by half in size, because half of the above area rectangle is empty.

The expected result would be George and and some "kid" George:



Yes I know, the syntax of the language can be improved but at least, everyone understands it. The most amusing in the above expression is that I could provide any picture to this higher order expression, and apply it to any rectangle, the behavior would be the same, the graphical result varying graphically speaking.

The only basic physical brick to my construction is the target rectangle. The following represents my rectangle abstraction:


A rectangle is made of three vectors. One representing the origin, the two others representing the horizontal and vertical directions.

Trust me when I tell you I tested the abstractions. You are going to see more tests. By the way in order not to clutter the article with too many pictures, I invite you to take some paper and reproduce the examples by hand, while reading the upcoming tests.

The abstractions involved into the rectangle reification process are

case class Rectangle (origin: Vector2D, horizontal: Vector2D, vertical: Vector2D)

case class Vector2D(x: Double, y: Double) {
  def +(other: Vector2D) = {
    Vector2D(x + other.x, y + other.y )
  }

  def * (scale: Double) = {
    Vector2D(scale * x, scale * y)
  }
}

I have extended my vector in order to add two vectors together and scale a vector when necessary. The tests are:

import org.specs2.mutable.Specification

final class TestVector2D  extends Specification{

  "A given vector " should {
    "bind its xcoordinate" in {
      Vector2D(3.0, 5.0).x must be equalTo (3.0)
    }

    "bind its ycoordinate" in {
      Vector2D(3.0, 5.0).x must be equalTo (3.0)
    }

  }

  "The addition of two vectors" should {
    " add the vectors coordinates " in {
      (Vector2D(3.0, 5.0) + Vector2D(7.0, 11.0)) must be equalTo(Vector2D(10.0, 16.0))
    }
  }

  "The scaling of a vector" should {
    " scale the  vector coordinates " in {
      (Vector2D(3.0, 5.0) * 2) must be equalTo(Vector2D(6.0, 10.0))
    }
  }
}

As I have planned to work with lines exclusively, I introduce a segment notion. A segment is defined by two vectors, one pointing to the start, the other pointing to the end:

case class Segment(start: Vector2D, end: Vector2D)

That's all. My incoming bricks are going to be functions and higher order functions. Recovering from many years of Java, I felt the need to introduce a "reference point" for the input picture to my chain of combinators. But what is a picture, if not a simple function aiming to take a rectangle as an input and draw something in it:

trait Picture extends Function1[Rectangle, Unit]

Dealing with lines, I defined a function LinearPicture, aiming to map a list of segments (a definition) into a target rectangle. Here comes the test:


import org.specs2.mutable.Specification

final class LinearPictureSpecification extends Specification{
    "Linear picture " should {
    "reproduce segment into unitary rectangle " in {
      implicit val drawer = new FakeDrawer
      val picture = LinearPicture(segments)(drawer)
      picture(inCell)
      drawer.segment must be equalTo (Segment(Vector2D(1, 1), Vector2D(2, 2)))
    }
  }

  def segments = List(Segment(Vector2D(0, 0), Vector2D(1, 1)))

  def inCell = Rectangle(Vector2D(1, 1), Vector2D(1, 0), Vector2D(0, 1))
}

The picture, will need on behalf the help of a drawer which knows how to draw physically lines, and a list of segments is provided as a definition for the graphic. Take five minutes with a paper and a pencil to reproduce the stuff. I will need a LineDrawer trait definition:

trait LineDrawer {
     def draw(segment: Segment)
}

and a little fake line drawer to trace what hes been drawn:

final class FakeDrawer extends LineDrawer{
  var segment: Segment = _

  def draw(aSegment: Segment) = segment = aSegment
}

Still with me? This is my implementation:

import com.promindis.geometry.FunctionalGeometry._

class LinearPicture(defintion: List[Segment], drawer: LineDrawer) extends Picture{

  def apply(rectangle: Rectangle) {
    val mapping = coordinateMappingIn(rectangle)
    defintion.foreach(segment =>
        drawer.draw(
          Segment(mapping(segment.start), mapping(segment.end))))
  }
}

object LinearPicture {
  def apply(defintion: List[Segment])(implicit drawer: LineDrawer): (Rectangle => Unit) = {
    new LinearPicture(defintion, drawer)
  }
}

The LinearPicture companion object can work with some implicit line drawer more suited to the context of execution. The most important job is achieved into the function apply method: iterating through each segment and applying the draw operation onto the rectangle. The apply method of the companion object, captures the definition of the graphic to be applied to whatever rectangle area we want

The mapping utility comes from the object FunctionalGeometry:

object FunctionalGeometry {

//..............

  def coordinateMappingIn(rectangle: Rectangle): (Vector2D => Vector2D) = {
    point: Vector2D =>
      rectangle.origin +
        rectangle.horizontal * point.x +
        rectangle.vertical * point.y
  }
//..............

We implicitly admitted that the genuine provided segment definitions will always be described by vectors scoped into the boundaries of a unitary rectangle.

From that point we have all we need. We can construct on our picture functionality. Take back your paper and pencil :)

You want to flip some picture function ? Then test the following:

import com.promindis.geometry.FunctionalGeometry._
import org.specs2.mutable.Specification

final class FlipSpecification extends Specification {

  "Flip a linear slope " should {
    "in a test rectangle mirror x coordinate versus vertical center axis" in {
      val drawer = new FakeDrawer
      val picture = LinearPicture(List(segment))(drawer)
      flip(picture)(inCell)
      drawer.segment must be equalTo (Segment(Vector2D(1.75, 1.25), Vector2D(1.25, 1.75) ))
    }
  }

  def segment = Segment(Vector2D(0.25, 0.25), Vector2D(0.75, 0.75))

  def inCell = Rectangle(Vector2D(1, 1), Vector2D(1, 0), Vector2D(0, 1))
}

and you get with:

object FunctionalGeometry {

  type G = Rectangle => Unit
//...........
  def flip(picture: G): G = {
    rectangle: Rectangle => {
      picture(Rectangle(
        rectangle.origin + rectangle.horizontal,
        rectangle.horizontal * -1,
        rectangle.vertical)
      )
    }
  }

//...........
}

Quite self explanatory. You want to set two pictures side by side with some ratio ?


Let's go for a test:

import com.promindis.geometry.FunctionalGeometry._
import org.specs2.mutable.Specification

final class BesideSpecification extends Specification {

  "Beside with two segments " should {
    "draw the two segemts in the box " in {
      implicit val drawer: DrawingAccumulator = new DrawingAccumulator()
      val p1 = LinearPicture(firstSegment)(drawer)
      val p2 = LinearPicture(secondSegment)(drawer)

      val sideBySide = beside(p1, p2, 0.5)
      sideBySide(inCell)
      println(drawer.segments)
      drawer.segments.contains(Segment(Vector2D(1, 1), Vector2D(1.5, 2))) must be equalTo (true)
      drawer.segments.contains(Segment(Vector2D(1.5, 2), Vector2D(2, 1))) must be equalTo (true)

    }
  }

  def firstSegment = List(Segment(Vector2D(0, 0), Vector2D(1, 1)))

  def secondSegment = List(Segment(Vector2D(0, 1), Vector2D(1, 0)))

  def inCell = Rectangle(Vector2D(1, 1), Vector2D(1, 0), Vector2D(0, 1))

}

Implementation follows immediately:

object FunctionalGeometry {

  type G = Rectangle => Unit
//...........

  def beside(pictureOne: G, pictureTwo: G, withRatio: Double): G = {
    rectangle: Rectangle => {
      pictureOne(inLeft(rectangle, withRatio))
      pictureTwo(inRight(rectangle, withRatio))
    }
  }

  private def inLeft(rectangle: Rectangle, ratio: Double): Rectangle = {
    Rectangle(rectangle.origin, rectangle.horizontal * ratio, rectangle.vertical)
  }

  private def inRight(rectangle: Rectangle, ratio: Double): Rectangle = {
    Rectangle(rectangle.origin + rectangle.horizontal * ratio, rectangle.horizontal * (1 - ratio), rectangle.vertical)
  }

}

And so on, we come with the small following language:

object FunctionalGeometry {

  type G = Rectangle => Unit


  def coordinateMappingIn(rectangle: Rectangle): (Vector2D => Vector2D) = {
    point: Vector2D =>
      rectangle.origin +
        rectangle.horizontal * point.x +
        rectangle.vertical * point.y
  }

  def flip(picture: G): G = {
    rectangle: Rectangle => {
      picture(Rectangle(
        rectangle.origin + rectangle.horizontal,
        rectangle.horizontal * -1,
        rectangle.vertical)
      )
    }
  }

  def empty(): G = {
    rectangle: Rectangle => ()
  }

  def mirror(picture: G): G = {
    rectangle: Rectangle => {
      picture(Rectangle(
        rectangle.origin + rectangle.vertical,
        rectangle.horizontal,
        rectangle.vertical * -1)
      )
    }
  }

  def rotate(picture: G): G = {
    rectangle: Rectangle => {
      picture(Rectangle(
        rectangle.origin + rectangle.vertical,
        rectangle.vertical * (-1),
        rectangle.horizontal)
      )
    }
  }

  def above(pictureOne: G, pictureTwo: G, withRatio: Double): G = {
    rectangle: Rectangle => {
      pictureOne(inTop(rectangle, withRatio))
      pictureTwo(inBottom(rectangle, withRatio))
    }
  }

  private def inTop(rectangle: Rectangle, ratio: Double): Rectangle = {
    Rectangle(rectangle.origin + rectangle.vertical * (1 - ratio), rectangle.horizontal, rectangle.vertical * (1 - ratio))
  }

  private def inBottom(rectangle: Rectangle, ratio: Double): Rectangle = {
    Rectangle(rectangle.origin, rectangle.horizontal, rectangle.vertical * ratio)
  }


  def beside(pictureOne: G, pictureTwo: G, withRatio: Double): G = {
    rectangle: Rectangle => {
      pictureOne(inLeft(rectangle, withRatio))
      pictureTwo(inRight(rectangle, withRatio))
    }
  }

  private def inLeft(rectangle: Rectangle, ratio: Double): Rectangle = {
    Rectangle(rectangle.origin, rectangle.horizontal * ratio, rectangle.vertical)
  }

  private def inRight(rectangle: Rectangle, ratio: Double): Rectangle = {
    Rectangle(rectangle.origin + rectangle.horizontal * ratio, rectangle.horizontal * (1 - ratio), rectangle.vertical)
  }

}

So what ? You want to see it ? Tough guys you are :). Ok, I took my Programming in Scala (2nd ed), and started to program a small GUI to display my pictures.

I came with a fast solution inspired by the book examples:

import scala.swing._
import java.awt.Graphics2D
import com.promindis.geometry.FunctionalGeometry._
import com.promindis.geometry.Schemes._
import com.promindis.graphics.GraphicsLineDrawer._
import com.promindis.geometry._
import com.promindis.graphics.Definitions._


class Painter(width: Int, height: Int, draw: G => G, draft: List[Segment]) extends Panel {
  def boundaries = Rectangle(Vector2D(0, 0), Vector2D(width.toDouble, 0), Vector2D(0, height.toDouble))

  preferredSize = new Dimension(width, height)
  border = Swing.EtchedBorder

  override def paintComponent(onGraphic: Graphics2D) {
    draw(LinearPicture(draft)(onGraphic))(boundaries)
  }
}

object DemoApp extends SimpleSwingApplication {
  def top = new MainFrame {
    title = "George and al."
    contents = new BoxPanel(Orientation.Vertical) {
      contents += new FlowPanel {
        contents += new Painter(500, 500, genuine(), george)
        contents += new Painter(500, 500, dancing, george)
      }
      contents += new FlowPanel {
        contents += new Painter(500, 500, bigGeorgeAndSmallGeorge(), george)
        contents += new Painter(500, 500, mosaic(), george)
      }
    }
    pack()
  }


  def genuine(): Scheme = {
    graphic: G => mirror(graphic)
  }

  def dancing = quartet()


  def bigGeorgeAndSmallGeorge() = {
    picture: G => mirror(beside(picture, above(empty(), picture, byHalf), byHalf))
  }

}

The class DemoApp extends SimpleSwingApplication and not SimpleGUIApplication as it has been deprecated since the book edition. I basically split the frame in four, playing with patterns. As I am not very found of Gui conception, I hope most of you will forgive my clumsy conception. The drawing job is achieve into the Painter class:

class Painter(width: Int, height: Int, draw: G => G, draft: List[Segment]) extends Panel {
  def boundaries = Rectangle(Vector2D(0, 0), Vector2D(width.toDouble, 0), Vector2D(0, height.toDouble))

  preferredSize = new Dimension(width, height)
  border = Swing.EtchedBorder

  override def paintComponent(onGraphic: Graphics2D) {
    draw(LinearPicture(draft)(onGraphic))(boundaries)
  }
}


where the boundaries of the rectangle area where to draw are set by the dimensions of the painter component . But by what magic, can the LinearPicture directly draw onto the Java2D graphics object? The mystery is solved when having a look at the following partner class definition:

import java.awt.{Color, Graphics2D}
class GraphicsLineDrawer(g: Graphics2D) extends LineDrawer{

  implicit def doubleToInt(value: Double): Int = value.toInt

  def draw(segment: Segment)  {
    g.setColor(Color.BLACK)
    g.drawLine(segment.start.x, segment.start.y, segment.end.x, segment.end.y)
  }
}

object GraphicsLineDrawer {

  implicit def graphics2DToGraphicsLineDrawer(g: Graphics2D): GraphicsLineDrawer = apply(g)

  def apply(g: Graphics2D) = {
    new GraphicsLineDrawer(g)
  }
}



I used the implicit conversion

implicit def graphics2DToGraphicsLineDrawer(g: Graphics2D): GraphicsLineDrawer = apply(g)

so to provide the illusion of a direct processing of the Java2D graphic object. No magic. The GraphicsLineDrawer "seems" to draw directly onto the Java2D graphic object instance:

draw(LinearPicture(draft)(onGraphic))(boundaries)

In order to play a little bit, you need to have some segment definitions:

object Definitions {
  val bigF = List(
    Segment(Vector2D(0.1, 0.1), Vector2D(0.3, 0.1)),
    Segment(Vector2D(0.3, 0.1), Vector2D(0.3, 0.3)),
    Segment(Vector2D(0.3, 0.3), Vector2D(0.7, 0.3)),
    Segment(Vector2D(0.7, 0.3), Vector2D(0.7, 0.5)),
    Segment(Vector2D(0.7, 0.5), Vector2D(0.3, 0.5)),
    Segment(Vector2D(0.3, 0.5), Vector2D(0.3, 0.7)),
    Segment(Vector2D(0.3, 0.7), Vector2D(0.9, 0.7)),
    Segment(Vector2D(0.9, 0.7), Vector2D(0.9, 0.9)),
    Segment(Vector2D(0.9, 0.9), Vector2D(0.1, 0.9)),
    Segment(Vector2D(0.1, 0.9), Vector2D(0.1, 0.1))
  )

  val george = List(
    Segment(Vector2D(0.3, 0), Vector2D(0.5, 0.4)),
    Segment(Vector2D(0.5, 0.4), Vector2D(0.6, 0)),
    Segment(Vector2D(0.8, 0), Vector2D(0.7, 0.5)),
    Segment(Vector2D(0.7, 0.5), Vector2D(1, 0.2)),
    Segment(Vector2D(1, 0.3), Vector2D(0.7, 0.7)),
    Segment(Vector2D(0.7, 0.7), Vector2D(0.5, 0.7)),
    Segment(Vector2D(0.5, 0.7), Vector2D(0.6, 0.8)),
    Segment(Vector2D(0.6, 0.8), Vector2D(0.5, 1)),
    Segment(Vector2D(0.4, 1), Vector2D(0.3, 0.8)),
    Segment(Vector2D(0.3, 0.8), Vector2D(0.4, 0.7)),
    Segment(Vector2D(0.4, 0.7), Vector2D(0.3, 0.7)),
    Segment(Vector2D(0.3, 0.7), Vector2D(0.2, 0.6)),
    Segment(Vector2D(0.2, 0.6), Vector2D(0, 0.8)),
    Segment(Vector2D(0, 0.7), Vector2D(0.1, 0.4)),
    Segment(Vector2D(0.1, 0.4), Vector2D(0.3, 0.6)),
    Segment(Vector2D(0.3, 0.6), Vector2D(0.4, 0.5)),
    Segment(Vector2D(0.4, 0.5), Vector2D(0.1, 0))
  )
}


and a complementary class, where I awkwardly started to defines "schemes" of behavior, as kind of patterns of graphics. This is where I will need some help, because I am still very low on combination and typing in Scala. But I am willing, so forgivable ;):

import com.promindis.geometry.FunctionalGeometry._

object Schemes {
  type Scheme = G => G

  def byHalf: Double = {
    0.5
  }

  def quartet (): Scheme = {
    picture: G => {
      val draw = pattern()
      above(draw(picture), mirror(draw(picture)), byHalf)
    }
  }

  def mosaic (): Scheme = {
    picture: G => {
      val draw = quartet()
      above(beside(draw(picture), draw(picture), byHalf), mirror(beside(draw(picture), draw(picture), byHalf)), byHalf)
    }
  }

  def pattern(): Scheme = (picture: G) => beside(picture, flip(picture), byHalf)

}

The running program provides the following picture:



and at a very low cost in time.

What we have built at this point is, layer after layer,

  • some bottom primitives including the notion of rectangles, segments, then pictures as functions, but seen as "data"
  • Then a geometric language above that layer of primitives
  • At the end we started a clumsy -but understandable-  language of schemes of combinations using the geometric layers

I can imagine that with some skilled Scala developers we could go further, and stack layers onto layers providing a nice complete embedded language.

Critics will be welcomed as this is an experimental work.
Time for Torchwood episode.


Be seeing you !!:)