[Help-gsl] Slow BLAS

[Gong Yi Liao]

Hello, 
   I find that the Linear algebra subroutines bundled (more specific,
the BLAS) with GSL1.8 are quite slow, for example, For computing SVD of
a matrix with high dimension (for example, 1000 by 1000), the GSL BLAS
(the official cblas, so called gslcblas) take about one minute (56.546
secs) on my laptop (Intel Duo Core T2300 with 768 ram and 768 swap
running Ubuntu Linux 6.06.1 LTS), mean while, within the popular
statistics programming environment, R, it just take less then ten
seconds (with GotoBLAS, ATLAS is slower than GotoBLAS without
significant difference) to compute SVD
the command:
> system.time(svd(matrix(rnorm(1e6), 1e3)))

the following is my code for computing SVD:
  /*
** gslbenc.c
** 
** Made by (Gong Yi Liao)
** Login   <address@hidden>
** 
** Started on  Sun Oct 22 13:14:41 2006 Gong Yi Liao
** Last update Sun Oct 22 13:56:52 2006 Gong Yi Liao
*/
#include <gsl/gsl_matrix.h>
#include <gsl/gsl_rng.h>
#include <gsl/gsl_randist.h>
#include <gsl/gsl_linalg.h>
/* #include <gsl/gsl_complex.h> */
/* #include <gsl/gsl_complex_math.h> */
/* #include <stdio.h> */
/* #include <stdlib.h> */


int main (void) {
  int i, j;
  const gsl_rng_type *T;
  gsl_rng *r;

  gsl_rng_env_setup();
  T = gsl_rng_default;
  r = gsl_rng_alloc (T);

  gsl_matrix *x = gsl_matrix_alloc (1000, 1000);
  
  // printf ("Start allocating matrix \n");
  for (i = 0; i < (x -> size1); i++) {
    for (j = 0; j < (x -> size2); j++)  {
      gsl_matrix_set(x, i, j, gsl_ran_gaussian(r, 1.0));
    }
  }
  //  printf ("matrix allocation finished\n"); 

  gsl_matrix *v = gsl_matrix_alloc (1000, 1000);
  gsl_vector *s = gsl_vector_alloc (1000);
  // gsl_matrix *mx = gsl_matrix_alloc (1000, 1000);
  gsl_vector *work = gsl_vector_alloc (1000);

  gsl_linalg_SV_decomp (x, v, s, work);

  // gsl_linalg_SV_decomp_mod  (x, mx, v, s, work);

  return 0;
}

with the gcc flags:
(with officical cblas)
gcc -o gslbenc -march=prescott -msse3 -mfpmath=sse,387 -I/usr/include
-lgsl -lm -lgslcblas gslbenc.c
(with ATLAS 3.7)
gcc -o gslbenc -march=prescott -msse3 -mfpmath=sse,387 -I/usr/include
-lgsl -lm -lcblas -latlas gslbenc.c
(with GotoBLAS)
gcc -o gslbenc -march=prescott -msse3 -mfpmath=sse,387 -I/usr/include
-lgsl -lm -lcblas -L/usr/lib/GotoBlas -lgoto gslbenc.c

all the binaries generated by these 3 compilation produce are slow SVD
executable, they take about 1 minute to compute SVD of a 1000 by 1000
matrix, and blas-based R just take 10 secs

it seems we need some huge refactoring the CBLAS part of GSL.