[Top][All Lists]
[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: [Help-gsl] Question about solving linear equation with GSL
From: |
Sheng Yu |
Subject: |
Re: [Help-gsl] Question about solving linear equation with GSL |
Date: |
Wed, 7 Jan 2009 14:55:36 -0500 |
On Wed, Jan 7, 2009 at 4:10 AM, Leo Razoumov <address@hidden> wrote:
> On 1/6/09, Sheng Yu <address@hidden> wrote:
>> Hi,
>>
>> I'm new to GSL. I'm trying to solve equation Ax=b, using one of the
>> solvers in Linear Algebra pack.
>> A is a 914x914 matrix. It is symmetric and positive definite. I use
>> Cholesky Decomposition to decomposite A first and then
>> solve it. Here is the pseudo code:
>>
>> const size_t N = 914;
>> gsl_matrix *A = gsl_matrix_alloc(N, N);
>> gsl_vector *b = gsl_vector_alloc(N);
>> gsl_vector *x = gsl_vector_alloc(N);
>>
>> // code to initialize A and b
>>
>> gsl_linalg_cholesky_decomp(A);
>> gsl_linalg_cholesky_solve(A, b, x);
>>
>> It takes 0.25 seconds for the code to do Cholesky Decomposition and
>> solve. However, when I do it with Matlab with A\b, it only
>> takes 0.045 seonds (It even does not take advantage of the fact that A
>> is symmetric and positive definite). I'm wondering why
>> Matlab is much faster. Did I do something wrong?
>>
>> Thank you!
>>
>> Yu Sheng
>>
>
> Few Questions:
> (1) What platform do you use (CPU type, Operation System, compiler, etc)?
I'm developing on an Intel core2 Quad Q9450 CPU. The operating system
is Fedora 9, and the compiler is gcc 4.3.0
> (2) What optimization flags do you use to compile GSL and your code.
> In the case of GCC "-O2" or "-O3" sould speed up things considerably.
> Dropping "-g" could also help.
I just turned on -O3 for optimization. But it seems that there is not
much difference with the running time no matter I turned it on or
not.
> (3) Is your matrix sparse?
No, it has a lot of 0s, but it's not sparse.
>
> --Leo--
>
Another library test I did is Intel Mathematics Kernel library. It
runs as fast as Matlab. The LAPACK library in it is implemented in
Fortran. I don't know whether it is exactly the same as the LAPACK
library privided by http://www.netlib.org/lapack/
Yu