Discriminated unions in Mercury

While reading the Mercury tutorial I learned how to define discriminated unions or algebraic datatypes.

Discriminated unions provide a way to define new structured data types.

For example here's a couple of type definitions that represent very simple arithmetic expressions:

:- type expression --->  const(int)
                   ; var(string)
                   ; binaryOperation(expression,operator,expression).
:- type operator ---> add ; subtract.

These definitions specify arithmetic expression parts such as constants(const), variables(var), binary operations(binaryOperation), etc.

Given these definitions we can instantiate them by calling the constructors:

main(!IO) :-
   X = binaryOperation(var("x"),add,const(10)),

Also we can use Mercury's unification to easily manipulate them, for example the following function writes the content of the expression to the output:

:- func operator_to_string(operator) = string.
operator_to_string(add)       = "+".
operator_to_string(subtract)  = "-".

:- pred write_expression(expression::in,io::di,io::uo) is det.

write_expression(const(C),!IO) :-
   io.write_int(C,!IO).
write_expression(var(S),!IO) :-
   io.write_string(S,!IO).
write_expression(binaryOperation(Operand1,Operator,Operand2) ,!IO) :-
   write_expression(Operand1,!IO),
   io.write_string(operator_to_string(Operator),!IO),
   write_expression(Operand2,!IO).

Calling this predicate like this:

main(!IO) :-
   X = binaryOperation(var("x"),add,const(10)),
   write_expression(X,!IO).

Generates:

x+10

Drawing the Apollonian Gasket with Common Lisp and Vecto

In this post I'm going to show a simple, yet incomplete, program to draw the Apollonian Gasket using Common Lisp and Vecto.

The Wikipedia entry for the Apollonian Gasket defines it as:

In mathematics, an Apollonian gasket or Apollonian net is a fractal generated from three circles, any two of which are tangent to one another. It is named after Greek mathematician Apollonius of Perga.

Basically the Apollonian gasket is generated by finding the 4th circle that is tangent to a given tree tangent circles and repeating the process recursively.

According to Wikipedia in order to find the position and the radius of the 4th circle we use the Descartes theorem which implies that we can apply the following formula:



Given that k_n is the curvature of the circle n. The curvature can be calculated as 1/radius of the circle.

According to the Wikipedia article this formula could also be used to determine the position of the 4th circle. By multiplying each curvature value k_n with a complex number that represents the position of the circle we can obtain a complex number with the position of the remaining circle.

For example given the following circles:



We can calculate the position of the 4th circle (shown in red) by using the following function:

(defun solve-equation (k1 k2 k3)
   (let ((equation 
           (lambda (op) 
                (+ k1 k2 k3 
                   (funcall op 
                       (* 2 (sqrt (+ (* k1 k2) 
                                     (* k2 k3) 
                                     (* k1 k3) ))))))))
       (values (funcall equation #'+)
               (funcall equation #'-))))

This function is used to solve the equation shown above. We use values to return multiple values implied by the solution of the equation . For this post I'm only using the first solution (that's why the program is incomplete!!!).

To find the radius we can say:

* (* 50 (solve-equation 1/50 1/50 1/50))

6.4641013

To find to position we use the complex numbers provided by Common Lisp and the same function defined above:

* (defvar radius  (* 50 (solve-equation 1/50 1/50 1/50)))

RADIUS
* (* radius  (solve-equation (* 1/50 #c(50 50)) (* 1/50 #c(150 50)) (* 1/50 #c(100 136.6))))

#C(83.56915 65.90832)

In order to draw the picture of the Apollonian gasket we use the Vecto library which provides an easy way for drawing the picture.

(defun draw-forth-circle (pos1 r1 pos2 r2 pos3 r3)
     (set-rgb-stroke  1 0 0)
     (let* ((curv1 (/ 1 r1))
            (curv2 (/ 1 r2))
            (curv3 (/ 1 r3))
            (kz1 (* curv1 pos1))
            (kz2 (* curv2 pos2))
            (kz3 (* curv3 pos3))
            (r4   (/ 1 (solve-equation curv1 curv2 curv3)))
            (pos4 (* r4 (solve-equation kz1 kz2 kz3))))
          (progn
             (centered-circle-path 
                    (realpart pos4) 
                    (imagpart pos4) 
                    r4)
             (stroke)
             (list pos4 r4))))


(defun draw-apollonian-gasket (c1 r1 c2 r2 c3 r3 steps)
    (if (> steps 0)
        (let ((lc41 (draw-forth-circle c1 r1 c2 r2 c3 r3)))
             (draw-apollonian-gasket c1 r1 c2 r2 (car lc41) (cadr lc41) (- steps 1))
             (draw-apollonian-gasket c2 r2 c3 r3 (car lc41) (cadr lc41) (- steps 1))
             (draw-apollonian-gasket c3 r3 c1 r1 (car lc41) (cadr lc41) (- steps 1)))))

(defun tria (file)
   (with-canvas (:width 300 :height 300)
     (set-rgb-stroke 1 1 1)
     (set-line-width 1)
     (translate 50  150) 

     (let* ((pos1 #c(50 50))
            (r1   50)
            (pos2 #c(150 50))
            (r2   50)
            (pos3 (complex 
                      100
                      (- 50 (sqrt (- (* 100 100) (* 50 50))))))
            (r3 50))

          (centered-circle-path (realpart pos1) (imagpart pos1) r1)
          (centered-circle-path (realpart pos2) (imagpart pos2) r2)
          (centered-circle-path (realpart pos3) (imagpart pos3) r3)
          (stroke)     
          (draw-apollonian-gasket 
               pos1 r1 
               pos2 r2
               pos3 r3 5))
     (save-png file)))

Running this program generates the following picture:



By applying scale and translate command we can obtain a better picture:



This program only draws one portion of the complete Apollonian Gasket, in order draw the complete picture, the other solutions of the equation shown above must be used. In future posts I'm going to try to complete the program.

The Mercury Programming Language

Mercury is a programming language that combines concepts from the logic and functional programming paradigms. In this post I'm going to show a couple of features of this language.

The hello world program looks like this:

:- module hello.


:- interface.
:- import_module io.
:- pred main(io::di, io::uo) is det.

:- implementation.

main(!IO) :-
   io.write_string("Hello world from Mercury!\n",!IO).

Where module hello is a module declaration, interface declares elements exported by this module and implementation contains the code for this module.

A couple of days ago I started reading the Ralph Becket's Mercury tutorial which is a very nice tutorial for learning the language. Also there's useful documentation available in the website including a Prolog to Mercury Transition Guide .

The language allows the creation of functions, for example here's the mandatory factorial function:

:- func factorial(int) = int.
factorial(N) = (if N = 0 then 1 else N*factorial(N-1)).

But also the language allows the creation of predicates that result in multiple solutions as in Prolog. For example see the following predicate that generates the numbers from 0 to N:

:- pred numbers(int::in, int::out) is nondet.

numbers(0,0).
numbers(N,R) :-
  N > 0,
    (
    numbers(N-1,R);
    R = N  
    ).

The pred section declares the numbers predicate. It says that the first argument is of type int and is for input and the second argument is also type int and is for output. Finally the is nondet says that this predicate may result in zero or more solutions.

The first rule says that when asking for numbers from 0, return 0. The second rule says that if asking for a number greater than 0 then return all the numbers from 0 to N-1 and then return N.

Using this predicate we can implement the following example from a previous post:

We need to write a program to find a 4-digit number such that ABCD*4 = DCBA where A,B,C,D are the digits of that number.

(Complete problem is described here)

A Mercury predicate that solves this problem looks like this:

:- pred solve_problem(list(int)::out) is nondet.

solve_problem(R) :-
  numbers(9,A), 
  numbers(9,B), 
  numbers(9,C), 
  numbers(9,D), 
  A > 0,
  (1000*A + 100*B + 10*C + D)*4 = (1000*D + 100*C + 10*B + A),
  R = [A,B,C,D].

As with Prolog the backtracking mechanism will find the values of A,B,C and D that satisfy the constraints.

In order to use the result of this predicate the solutions predicate is used. This predicate returns a list with all the results from the specified predicate.

main(!IO) :-
      solutions(solve_problem,R),
      io.write(R,!IO),
      io.write_string("\n",!IO).

This program writes:


[[2, 1, 7, 8]]


In future post I'll try to explore more features of this interesting language.

Starting with Vecto

A couple of days ago I read about a library for Common Lisp called Vecto that allows the creation of graphics by describing them using Lisp functions.

Before trying this library, I noticed that it could be installed by using ASDF-Install. ASDF is like a packaging system for Common Lisp libraries and ASDF-Install is a program that let's you locate and install thouse libraries with its dependencies.

I'm using Steel Bank Common Lisp which already has ASDF and ASDF-Install integrated with it.

In order to install Vecto the only thing that I did was:

* (require :asdf)

NIL
* (require :asdf-install)
...
* (asdf-install:install :Vecto)

After answering a couple of questions It downloaded and installed Vecto with all its dependencies.


The library has nice documentation. Here's a little example I created for this post.

(require :vecto)
(defpackage #:vecto-test
  (:use #:cl #:vecto))

(in-package #:vecto-test)

(defun triangle (s c)
    (set-rgb-stroke (car c)  
                    (cadr c) 
                    (caddr c))
    (set-line-width 10)
    (set-line-join :round)
    (move-to 0 0)
    (line-to (- s) 0)
    (line-to 0 (- s))
    (line-to s 0)
    (close-subpath)
    (stroke)
 )

(defun gen-values (from to step)
  (if (<= from to)
     (cons from 
           (gen-values (+ from step) to step))
     nil))

(defun tria (file)
   (with-canvas (:width 300 :height 300)
     (mapcar #'(lambda (x)
      (with-graphics-state
         (translate 150 150)
         (rotate (cadr x))
         (triangle (car x) (list (random 1.0) 
                            (random 1.0) 
                            (random 1.0)))))
      (mapcar #'(lambda (length angle) (list length angle))
              (gen-values 10 200 30)
              (gen-values 0.0 3.14 0.4)))        
     (save-png file)))

(tria "vecto-test-image.png")

This code generates the following image: