gee-guile-gsl


GEE/Guile-GSL tutorial

Guile-GSL is a C/Scheme library extension for Guile, the GNU's Ubiquitous Intelligent Language for Extensions. It implements a binding to the GNU Scientific Library (GSL). The purpose is to implement an environment similar to GNU Octave, but using Scheme as a programming language.

Guile-GSL is a part of GEE, a collection of extensions for Guile.

This is a tutorial on the usage of Guile-GSL, not on how to program with Scheme. To understand what's going on here you have to know the basics of Scheme: how to define variables, how to define procedures, how to invoke procedures, what is a cons.

Many people can learn the Scheme basics in a few days; the Internet is full of introduction papers on the subject. The first chapters of the SICP are perfect:

 Structure and Interpretation of Computer Programs (sicp)
 Harold Abelson and Gerald Sussman
 Copyright (C) 1996 Massachussets Institute of Technology

this book is available on the Internet, there is also a GNU Texinfo version.

On this wiki are available installation instructions GEE-Guile-GSL-install.

Loading the module

To load the module:

 (define-module (this-script) 
   #:use-module (oop goops) 
   #:use-module (gee math oop) 
   #:use-module (gee math gsl) 
   #:duplicates merge-generics) 
  
 ;; your code goes here 

Interactive shell

When Guile-GSL is installed a couple of scripts are placed on the system to use the modules from the Guile REPL (read, eval, print loop):

gsh

a Bourne shell interface scripts that sets up the environment and then runs guile with gsh.scm as script; installed under a path like /usr/local/bin/gsh;

gsh.scm

a Guile script that loads Guile-GSL and other modules then invokes the Guile REPL; it is installed under a path like /usr/local/libexec/guile-gsl/0.1.0/gsh.scm.

If the installation was successful at the shell prompt we can do:

$ gsh
Guile-GSL version 0.1.0
Written by Marco Maggi.
     
Copyright (C) 2006, 2007 by Marco Maggi.
     
This is free software; see the  source or use the '--license' option for
copying conditions.  There is NO warranty; not  even for MERCHANTABILITY
or FITNESS FOR A PARTICULAR PURPOSE.
     
gsh>

gsh> is the default Guile-GSL shell prompt.

Both the scripts are very simple and you are encouraged to inspect and modify them to suit your needs. The gsh shell script is there exactly to be customised with site-specific environment variables.

Vectors and matrices

Vectors and matrices are the basic objects of math computations with Guile-GSL. The easiest interface is the one that uses the GOOPS extension of Guile; there is only one class that we have to be aware of: <gsl>, and we rarely need to use it explicitly.

Quick code to overview:

 (define a #,(gvr 3  1 2 3)) 
 (define b #,(gvr 3  4 5 6)) 
  
       
 (display (+ a b)) ;; -> [5 7 9] 
 (write   (+ a b)) ;; -> #,(gvr 3  5 7 9) 
       
 (display (dot   a b)) ;; -> 32.0 
 (display (cross a b)) ;; -> [-3.0 6.0 -3.0] 

We may want to define vectors and matrices in one of three ways: just define an "empty" object, define an object with known constant elements, define an object with elements that we have to compute using math functions.

Constant elements

We use the Scheme "read syntax" with the SRFI-10 features included in Guile. There are 4 symbols:

gvr

to declare vectors of real numbers (gvr = GSL vector real);

gvc

to declare vectors of complex numbers (gvc = GSL vector complex);

gmr

to declare matrices of real numbers (gmr = GSL matrix real);

gmc

to declare matrices of complex numbers (gmc = GSL matrix complex).

the definitions look like the following:

 #,(gvr 3  1 2 3) 
       
 #,(gvc 3  1+2i 2-3i 3+4i) 
       
 #,(gmr (2 . 3)  1 2 3 
                  4 5 6) 
       
 #,(gmc (2 . 3)  1+2i 2-3i 3+4i 
                  4-5i 6+7i 8-9i) 

The first element after the tag must be the dimension: for vectors the number of elements; for matrices a cons holding the number of rows in the car and the number of columns in the cdr. Notice that when using the read syntax the cons does not need to be quoted.

After the tag and the dimension we simply write the elements separated by blanks (spaces, tabs, newlines). The elements must be numbers.

Empty objects

Sometimes we need to define an object with known dimension that will be filled with values later. We do it like this:

 (define a (make-gsl-vector-real    3)) 
 (define b (make-gsl-vector-complex 3)) 
 (define c (make-gsl-matrix-real    '(3 . 5))) 
 (define d (make-gsl-matrix-complex '(19 . 2))) 

as Scheme dictates: outside the read syntax, the conses must be quoted. All the elements are set to 0.0; optionally we can give an additional parameter representing the default:

 (define a (make-gsl-vector-real 3  4.5)) 
 ;; -> [4.5 4.5 4.5] 

A special maker exists for identity matrices:

 (eye 3) 
 ;; -> [[1 0 0] 
 ;;     [0 1 0] 
 ;;     [0 0 1]] 

Objects with computed values

Defining objects with computed values is similar to defining empty objects in that we invoke functions for it:

 (define a (read-gsl-vector-real 3 
             (sin 0.2) (cos 0.5) (tan 0.6))) 
       
 (define b (read-gsl-vector-complex 3 
             (sin 0.2+9i) (cos 0.5) (tan -0.6i))) 
       
 (define c (read-gsl-matrix-real '(2 . 3) 
             (sin 0.2) (cos 0.5) (tan 0.6) 
             (sin 0.3) (cos 0.9) (tan 0.2))) 
       
 (define d (read-gsl-matrix-complex '(2 . 3) 
             (sin 0.2)    0.5       (tan 0.6) 
             (sin 0.2+9i) (cos 0.5) (tan -0.6i))) 

Functions

Once we have defined vectors and matrices we just apply functions to them. The four arithmetic functions +, -, * and / do what we expect element by element; there is also the unary - that negates all the elements:

 (- #,(gvr 3  1 2 3)) 
 ;; -> [-1 -2 -3] 

Other functions, like hyperbolic ones, are applied to all the elements, the following form:

 (sinh #,(gvr 3  1 2 3)) 

is equivalent to:

 (read-gsl-vector-real 3  (sinh 1) (sinh 2) (sinh 3)) 

DOT does the row by column product; if the arguments are both vectors: the operation is the scalar product; if the arguments are a matrix and a vector:

 (dot #,(gmc (2 . 3)  1+2i 2-3i 3+4i 
                      4-5i 6+7i 8-9i) 
      #,(gvr 3  1 2 3)) 
 ;; -> #,(gvc 2 14.0+8.0i 40.0-18.0i) 

and if they are two matrices:

 (dot #,(gmc (2 . 3)  1+2i 2-3i 3+4i 
                      4-5i 6+7i 8-9i) 
      #,(gmc (3 . 2)  1+2i 2-3i 
                      3+4i 4-5i 
                      6+7i 8-9i)) 
 ;; -> #,(gmc (2 . 2) 5.0+48.0i   61.0-16.0i 
 ;;                   115.0+50.0i 35.0-168.0i) 

There is no such concept as row vector or column vector: the vectors are interpreted as rows or columns as need be.

We can transpose matrices like this:

 (transpose m) 

and do the conjugate:

 (conjugate m) 

Setters and getters

The setter/getter synopsis is ugly, but that's the way it is. The synopsis of all the getters is:

 (type object key) 

so to get elements from a vector:

 (define a #,(gvr 3  1 2 3)) 
 (elm a 0) ;; -> 1 
 (elm a 1) ;; -> 2 
  
 (elm a 2) ;; -> 3 

indexes are zero based. To get elements from a matrix we have to select both the row and the column:

 (define a #,(gmr (2 . 3)  1 2 3 
                           4 5 6)) 
 (elm a '(0 . 0)) ;; -> 1 
 (elm a '(0 . 1)) ;; -> 2 
 (elm a '(1 . 2)) ;; -> 6 

to get a whole row:

 (row a 0) ;; -> [1 2 3] 
 (row a 1) ;; -> [4 5 6] 

and to get a whole column:

 (column a 0) ;; -> [1 4] 
 (column a 1) ;; -> [2 5] 

The synopsis of all the setters is:

 (set! (type object key) value) 

to set elements:

 (define a #,(gvr 3  1 2 3)) 
 (set! (elm a 0) -1) ;; -> [-1 2 3] 
       
 (define b #,(gmr (2 . 3)  1 2 3 
                           4 5 6)) 
 (set! (elm b '(1 . 2)) -1) 
 ;; -> #,(gmr (2 . 3)  1 2 3 
 ;;                    4 5 -1) 

to set a whole row or column:

 (define b #,(gmr (2 . 3)  1 2 3 
                           4 5 6)) 
       
 (set! (row b 0) #,(gvr 3  -1 -2 -3)) 
 ;; -> #,(gmr (2 . 3)  -1 -2 -3 
 ;;                    4 5 6) 
       
 (set! (column b 0) #,(gvr 3  -8 -9)) 
 ;; -> #,(gmr (2 . 3)  -8 -2 -3 
 ;;                    -9 5 6) 

to set a diagonal:

 (define b #,(gmr (2 . 3)  1 2 3 
                           4 5 6)) 
       
 (set! (diagonal b 1) #,(gvr 2  -8 -9)) 
 ;; -> #,(gmr (2 . 3)  1 -8 3 
 ;;                    4 5 -9) 

category-gee