%perlcode %{

use strict;
use warnings;
use Carp qw/croak/;
use Scalar::Util 'blessed';

use Math::GSL qw/:all/;
use Math::GSL::Errno    qw/:all/;
use Math::GSL::Eigen    qw/:all/;
use Math::GSL::Vector   qw/:all/;
use Math::GSL::Complex  qw/:all/;
use Math::GSL::Test  qw/is_similar/;
use Math::GSL::BLAS  qw/gsl_blas_dgemm/;
use Math::GSL::CBLAS qw/$CblasNoTrans/;
use Math::GSL::Permutation;
use Math::GSL::Linalg qw/
    gsl_linalg_LU_decomp
    gsl_linalg_LU_det
    gsl_linalg_LU_lndet
    gsl_linalg_LU_invert
/;

# should only include needed methods
use Math::GSL::MatrixComplex qw/:all/;
use Math::GSL::VectorComplex qw/:all/;
use Data::Dumper;

use overload
    '*'      => \&_multiplication,
    '+'      => \&_addition,
    '-'      => \&_subtract,
    '=='     => \&_equal,
    '!='     => \&_not_equal,
    fallback => 1;

our @EXPORT_OK = qw/
        gsl_matrix_alloc gsl_matrix_calloc gsl_matrix_alloc_from_block
        gsl_matrix_alloc_from_matrix gsl_vector_alloc_row_from_matrix
        gsl_vector_alloc_col_from_matrix gsl_matrix_free gsl_matrix_submatrix
        gsl_matrix_row gsl_matrix_column gsl_matrix_diagonal
        gsl_matrix_subdiagonal gsl_matrix_superdiagonal gsl_matrix_subrow
        gsl_matrix_subcolumn gsl_matrix_view_array
        gsl_matrix_view_array_with_tda gsl_matrix_view_vector
        gsl_matrix_view_vector_with_tda gsl_matrix_const_submatrix
        gsl_matrix_const_row gsl_matrix_const_column gsl_matrix_const_diagonal
        gsl_matrix_const_subdiagonal gsl_matrix_const_superdiagonal
        gsl_matrix_const_subrow gsl_matrix_const_subcolumn
        gsl_matrix_const_view_array gsl_matrix_const_view_array_with_tda
        gsl_matrix_const_view_vector gsl_matrix_const_view_vector_with_tda
        gsl_matrix_get gsl_matrix_set gsl_matrix_ptr gsl_matrix_const_ptr
        gsl_matrix_set_zero gsl_matrix_set_identity gsl_matrix_set_all
        gsl_matrix_fread gsl_matrix_fwrite gsl_matrix_fscanf gsl_matrix_fprintf
        gsl_matrix_memcpy gsl_matrix_swap gsl_matrix_swap_rows
        gsl_matrix_swap_columns gsl_matrix_swap_rowcol gsl_matrix_transpose
        gsl_matrix_transpose_memcpy gsl_matrix_max gsl_matrix_min
        gsl_matrix_minmax gsl_matrix_max_index gsl_matrix_min_index
        gsl_matrix_minmax_index gsl_matrix_isnull gsl_matrix_ispos
        gsl_matrix_isneg gsl_matrix_isnonneg gsl_matrix_add gsl_matrix_sub
        gsl_matrix_mul_elements gsl_matrix_div_elements gsl_matrix_scale
        gsl_matrix_add_constant gsl_matrix_add_diagonal
        gsl_matrix_char_alloc gsl_matrix_char_calloc  gsl_matrix_char_alloc_from_block
        gsl_matrix_char_alloc_from_matrix gsl_vector_char_alloc_row_from_matrix gsl_vector_char_alloc_col_from_matrix
        gsl_matrix_char_free gsl_matrix_char_submatrix
        gsl_matrix_char_row  gsl_matrix_char_column
        gsl_matrix_char_diagonal gsl_matrix_char_subdiagonal gsl_matrix_char_superdiagonal
        gsl_matrix_char_subrow gsl_matrix_char_subcolumn gsl_matrix_char_view_array
        gsl_matrix_char_view_array_with_tda gsl_matrix_char_view_vector gsl_matrix_char_view_vector_with_tda
        gsl_matrix_char_const_submatrix gsl_matrix_char_const_row gsl_matrix_char_const_column
        gsl_matrix_char_const_diagonal gsl_matrix_char_const_subdiagonal gsl_matrix_char_const_superdiagonal
        gsl_matrix_char_const_subrow gsl_matrix_char_const_subcolumn gsl_matrix_char_const_view_array
        gsl_matrix_char_const_view_array_with_tda gsl_matrix_char_const_view_vector gsl_matrix_char_const_view_vector_with_tda
        gsl_matrix_char_get gsl_matrix_char_set gsl_matrix_char_ptr gsl_matrix_char_const_ptr
        gsl_matrix_char_set_zero gsl_matrix_char_set_identity
        gsl_matrix_char_set_all  gsl_matrix_char_fread
        gsl_matrix_char_fwrite gsl_matrix_char_fscanf gsl_matrix_char_fprintf
        gsl_matrix_char_memcpy gsl_matrix_char_swap
        gsl_matrix_char_swap_rows gsl_matrix_char_swap_columns
        gsl_matrix_char_swap_rowcol gsl_matrix_char_transpose gsl_matrix_char_transpose_memcpy
        gsl_matrix_char_max gsl_matrix_char_min
        gsl_matrix_char_minmax  gsl_matrix_char_max_index
        gsl_matrix_char_min_index gsl_matrix_char_minmax_index
        gsl_matrix_char_isnull gsl_matrix_char_ispos gsl_matrix_char_isneg
        gsl_matrix_char_isnonneg  gsl_matrix_char_add
        gsl_matrix_char_sub gsl_matrix_char_mul_elements gsl_matrix_char_div_elements
        gsl_matrix_char_scale gsl_matrix_char_add_constant gsl_matrix_char_add_diagonal
        gsl_matrix_int_alloc gsl_matrix_int_calloc gsl_matrix_int_alloc_from_block
        gsl_matrix_int_alloc_from_matrix gsl_vector_int_alloc_row_from_matrix gsl_vector_int_alloc_col_from_matrix
        gsl_matrix_int_free gsl_matrix_int_submatrix gsl_matrix_int_row
        gsl_matrix_int_column gsl_matrix_int_diagonal gsl_matrix_int_subdiagonal
        gsl_matrix_int_superdiagonal gsl_matrix_int_subrow gsl_matrix_int_subcolumn gsl_matrix_int_view_array
        gsl_matrix_int_view_array_with_tda gsl_matrix_int_view_vector gsl_matrix_int_view_vector_with_tda
        gsl_matrix_int_const_submatrix gsl_matrix_int_const_row gsl_matrix_int_const_column
        gsl_matrix_int_const_diagonal gsl_matrix_int_const_subdiagonal gsl_matrix_int_const_superdiagonal
        gsl_matrix_int_const_subrow gsl_matrix_int_const_subcolumn gsl_matrix_int_const_view_array
        gsl_matrix_int_const_view_array_with_tda gsl_matrix_int_const_view_vector gsl_matrix_int_const_view_vector_with_tda
        gsl_matrix_int_get gsl_matrix_int_set
        gsl_matrix_int_ptr gsl_matrix_int_const_ptr
        gsl_matrix_int_set_zero gsl_matrix_int_set_identity gsl_matrix_int_set_all
        gsl_matrix_int_fread gsl_matrix_int_fwrite
        gsl_matrix_int_fscanf gsl_matrix_int_fprintf
        gsl_matrix_int_memcpy gsl_matrix_int_swap
        gsl_matrix_int_swap_rows gsl_matrix_int_swap_columns gsl_matrix_int_swap_rowcol
        gsl_matrix_int_transpose gsl_matrix_int_transpose_memcpy
        gsl_matrix_int_max gsl_matrix_int_min gsl_matrix_int_minmax
        gsl_matrix_int_max_index gsl_matrix_int_min_index
        gsl_matrix_int_minmax_index  gsl_matrix_int_isnull
        gsl_matrix_int_ispos gsl_matrix_int_isneg gsl_matrix_int_isnonneg
        gsl_matrix_int_add gsl_matrix_int_sub
        gsl_matrix_int_mul_elements gsl_matrix_int_div_elements gsl_matrix_int_scale
        gsl_matrix_int_add_constant gsl_matrix_int_add_diagonal
        gsl_matrix_get_row gsl_matrix_get_col gsl_matrix_set_row gsl_matrix_set_col
        /;

our %EXPORT_TAGS = ( all => [ @EXPORT_OK ],
                 char => [ qw/
                    gsl_matrix_char_alloc
                    gsl_matrix_char_calloc
                    gsl_matrix_char_alloc_from_block
                    gsl_matrix_char_alloc_from_matrix
                    gsl_vector_char_alloc_row_from_matrix
                    gsl_vector_char_alloc_col_from_matrix
                    gsl_matrix_char_free
                    gsl_matrix_char_submatrix
                    gsl_matrix_char_row
                    gsl_matrix_char_column
                    gsl_matrix_char_diagonal
                    gsl_matrix_char_subdiagonal
                    gsl_matrix_char_superdiagonal
                    gsl_matrix_char_subrow
                    gsl_matrix_char_subcolumn
                    gsl_matrix_char_view_array
                    gsl_matrix_char_view_array_with_tda
                    gsl_matrix_char_view_vector
                    gsl_matrix_char_view_vector_with_tda
                    gsl_matrix_char_const_submatrix
                    gsl_matrix_char_const_row
                    gsl_matrix_char_const_column
                    gsl_matrix_char_const_diagonal
                    gsl_matrix_char_const_subdiagonal
                    gsl_matrix_char_const_superdiagonal
                    gsl_matrix_char_const_subrow
                    gsl_matrix_char_const_subcolumn
                    gsl_matrix_char_const_view_array
                    gsl_matrix_char_const_view_array_with_tda
                    gsl_matrix_char_const_view_vector
                    gsl_matrix_char_const_view_vector_with_tda
                    gsl_matrix_char_get
                    gsl_matrix_char_set
                    gsl_matrix_char_ptr
                    gsl_matrix_char_const_ptr
                    gsl_matrix_char_set_zero
                    gsl_matrix_char_set_identity
                    gsl_matrix_char_set_all
                    gsl_matrix_char_fread
                    gsl_matrix_char_fwrite
                    gsl_matrix_char_fscanf
                    gsl_matrix_char_fprintf
                    gsl_matrix_char_memcpy
                    gsl_matrix_char_swap
                    gsl_matrix_char_swap_rows
                    gsl_matrix_char_swap_columns
                    gsl_matrix_char_swap_rowcol
                    gsl_matrix_char_transpose
                    gsl_matrix_char_transpose_memcpy
                    gsl_matrix_char_max
                    gsl_matrix_char_min
                    gsl_matrix_char_minmax
                    gsl_matrix_char_max_index
                    gsl_matrix_char_min_index
                    gsl_matrix_char_minmax_index
                    gsl_matrix_char_isnull
                    gsl_matrix_char_ispos
                    gsl_matrix_char_isneg
                    gsl_matrix_char_isnonneg
                    gsl_matrix_char_add
                    gsl_matrix_char_sub
                    gsl_matrix_char_mul_elements
                    gsl_matrix_char_div_elements
                    gsl_matrix_char_scale
                    gsl_matrix_char_add_constant
                    gsl_matrix_char_add_diagonal
                     /],

                 double => [ qw/
                    gsl_matrix_alloc
                    gsl_matrix_calloc
                    gsl_matrix_alloc_from_block
                    gsl_matrix_alloc_from_matrix
                    gsl_vector_alloc_row_from_matrix
                    gsl_vector_alloc_col_from_matrix
                    gsl_matrix_free
                    gsl_matrix_submatrix
                    gsl_matrix_row
                    gsl_matrix_column
                    gsl_matrix_diagonal
                    gsl_matrix_subdiagonal
                    gsl_matrix_superdiagonal
                    gsl_matrix_subrow
                    gsl_matrix_subcolumn
                    gsl_matrix_view_array
                    gsl_matrix_view_array_with_tda
                    gsl_matrix_view_vector
                    gsl_matrix_view_vector_with_tda
                    gsl_matrix_const_submatrix
                    gsl_matrix_const_row
                    gsl_matrix_const_column
                    gsl_matrix_const_diagonal
                    gsl_matrix_const_subdiagonal
                    gsl_matrix_const_superdiagonal
                    gsl_matrix_const_subrow
                    gsl_matrix_const_subcolumn
                    gsl_matrix_const_view_array
                    gsl_matrix_const_view_array_with_tda
                    gsl_matrix_const_view_vector
                    gsl_matrix_const_view_vector_with_tda
                    gsl_matrix_get
                    gsl_matrix_set
                    gsl_matrix_ptr
                    gsl_matrix_const_ptr
                    gsl_matrix_set_zero
                    gsl_matrix_set_identity
                    gsl_matrix_set_all
                    gsl_matrix_fread
                    gsl_matrix_fwrite
                    gsl_matrix_fscanf
                    gsl_matrix_fprintf
                    gsl_matrix_memcpy
                    gsl_matrix_swap
                    gsl_matrix_swap_rows
                    gsl_matrix_swap_columns
                    gsl_matrix_swap_rowcol
                    gsl_matrix_transpose
                    gsl_matrix_transpose_memcpy
                    gsl_matrix_max
                    gsl_matrix_minmax
                    gsl_matrix_max_index
                    gsl_matrix_min_index
                    gsl_matrix_minmax_index
                    gsl_matrix_isnull
                    gsl_matrix_ispos
                    gsl_matrix_isneg
                    gsl_matrix_isnonneg
                    gsl_matrix_add
                    gsl_matrix_mul_elements
                    gsl_matrix_div_elements
                    gsl_matrix_scale
                    gsl_matrix_add_constant
                    gsl_matrix_add_diagonal
                    /],
int => [ qw/
                 gsl_matrix_int_alloc
                 gsl_matrix_int_alloc_from_matrix
                 gsl_matrix_int_free
                 gsl_matrix_int_column
                 gsl_matrix_int_superdiagonal
                 gsl_matrix_int_view_array_with_tda
                 gsl_matrix_int_const_submatrix
                 gsl_matrix_int_const_diagonal
                 gsl_matrix_int_const_subrow
                 gsl_matrix_int_const_view_array_with_tda
                 gsl_matrix_int_get
                 gsl_matrix_int_ptr
                 gsl_matrix_int_set_zero
                 gsl_matrix_int_fread
                 gsl_matrix_int_fscanf
                 gsl_matrix_int_memcpy
                 gsl_matrix_int_swap_rows
                 gsl_matrix_int_transpose
                 gsl_matrix_int_max
                 gsl_matrix_int_max_index
                 gsl_matrix_int_minmax_index
                 gsl_matrix_int_ispos
                 gsl_matrix_int_add
                 gsl_matrix_int_mul_elements
                 gsl_matrix_int_add_constant
                 /],
);


=encoding utf8

=head1 NAME

Math::GSL::Matrix - Mathematical functions concerning Matrices

=head1 SYNOPSIS

    use Math::GSL::Matrix qw/:all/;
    my $matrix1 = Math::GSL::Matrix->new(5,5);  # OO interface
    my $matrix2 = $matrix1 + 4;                 # You can add or subtract values or matrices to OO matrices
    my $matrix3 = $matrix1 - 4;
    my $matrix4 = $matrix2 + $matrix1;
    my $matrix5 = $matrix2 . $matrix1;          # This is a scalar product, it simply multiply each element
                                                # with the element of $matrix1 that have the same position
                                                # See Math::GSL::BLAS if you want scalar product

    my $matrix6 = $matrix2 . 8;                 # Multiply every elements of $matrix2 by 8
    my $matrix7 = $matrix2 * $matrix1;          # scalar product of two matrices
    if($matrix1 == $matrix4) ...
    if($matrix1 != $matrix3) ...
    my $matrix8 = gsl_matrix_alloc(5,5);        # standard interface

=head1 DESCRIPTION

This module is part of the L<Math::GSL> distribution. It defines a Perl insterface
to GNU Scientific Library matrices.

There are two different (but not exclusive) ways to use this module: using the
OO API, built  manually over the GSL functions, or using directly the
functions defined by GSL library.

=head1 OBJECT ORIENTED API

=head2 Constructor

=head3 Math::GSL::Matrix->new()

Creates a new Matrix of the given size.

    my $matrix = Math::GSL::Matrix->new(10,10);

If by any chance you already have a gsl_matrix and want to "objectify" it, you
can use the following constructor:

    my $m = gsl_matrix_alloc(10, 20);
    # ... something
    my $matrix = Math::GSL::Matrix->new($m);

Note that C<$m> is B<NOT> copied. The new object will just refer to the old gsl_matrix.

=cut

sub new {
    return _new_from_matrix(@_) if @_ == 2;
    return _new_from_dims(@_)   if @_ == 3;

    croak( __PACKAGE__ .
          sprintf("::new(...) - unknown Matrix constructor with %d args", scalar(@_)));
}

sub _new_from_matrix
{
    my ($class, $m) = @_;
    if (!defined($m) || !blessed($m) || ref($m) ne "Math::GSL::Matrix::gsl_matrix") {
        croak( __PACKAGE__.'::new($m) - $m should be a gsl_matrix');
    }

    return bless {
        _rows => $m->swig_size1_get(),
        _cols => $m->swig_size2_get(),
        _matrix => $m,
    }, $class;
}

sub _new_from_dims
{
    my ($class, $rows, $cols) = @_;
    my $this = {};
    my $matrix;
    if ( defined $rows       && defined $cols &&
        $rows > 0            && $cols > 0     &&
        (int $rows == $rows) && (int $cols == $cols)){

        $matrix  = gsl_matrix_alloc($rows,$cols);
    } else {
        croak( __PACKAGE__.'::new($x,$y) - $x and $y must be positive integers');
    }
    gsl_matrix_set_zero($matrix);
    $this->{_matrix} = $matrix;
    ($this->{_rows}, $this->{_cols}) = ($rows,$cols);
    bless $this, $class;
}

=head2 Getters

=head3 raw()

Get the underlying GSL matrix object created by SWIG, useful for using gsl_matrix_* functions which do not have an OO counterpart.

    my $matrix     = Math::GSL::Matrix->new(3,3);
    my $gsl_matrix = $matrix->raw;
    my $stuff      = gsl_matrix_get($gsl_matrix, 1, 2);

=cut

sub raw  { (shift)->{_matrix} }


=head3 dim()

Returns the number of rows and columns in a matrix as an array.

   my ($rows, $cols) = $matrix->dim;

Basically a shortcut to C<rows> and C<cols> methods.

=cut

sub dim { ($_[0]->rows, $_[0]->cols) }

=head3 rows()

Returns the number of rows in the matrix.

    my $rows = $matrix->rows;
=cut

sub rows { (shift)->{_rows}   }

=head3 cols()

Returns the number of columns in the matrix.

    my $cols = $matrix->cols;
=cut

sub cols { (shift)->{_cols}   }

=head3 as_list()

Get the contents of a Math::GSL::Matrix object as a Perl list.

    my $matrix = Math::GSL::Matrix->new(3,3);
    ...
    my @matrix = $matrix->as_list;
=cut

sub as_list($)
{
    my $self = shift;
    my @matrix;
    for my $row ( 0 .. $self->rows-1) {
       push @matrix, map { gsl_matrix_get($self->raw, $row, $_) } (0 .. $self->cols-1 );
    }
    return @matrix;
}

=head3 as_vector()

Returns a 1xN or Nx1 matrix as a Math::GSL::Vector object. Dies if called on a matrix that is not a single row or column. Useful for turning the output of C<col()> or C<row()> into a vector, like so:

    my $vector1 = $matrix->col(0)->as_vector;
    my $vector2 = $matrix->row(1)->as_vector;

=cut

sub as_vector($)
{
    my ($self) = @_;
    croak(__PACKAGE__.'::as_vector() - must be a single row or column matrix') unless ($self->cols == 1 or $self->rows == 1);

    # TODO: there is a faster way to do this
    return Math::GSL::Vector->new([ $self->as_list ] );

}

=head3 get_elem()

Returns an element of a matrix.

    my $matrix = Math::GSL::Matrix->new(3,3);
    ...
    my $middle = $matrix->get_elem(1,1);

B<NOTE:> just like any other method on this module, rows and arrays start
with indice 0.

=cut

sub get_elem {
    my ($self, $row, $col) = @_;
    die __PACKAGE__.'::get_elem($x, $y, $v) - $x must be a valid row number'
        if ($row < 0 || $row >= $self->rows);
    die __PACKAGE__.'::get_elem($x, $y, $v) - $y must be a valid column number'
        if ($col < 0 || $col >= $self->cols);
    return gsl_matrix_get($self->raw, $row, $col);
}

=head3 row()

Returns a row matrix of the row you enter.

    my $matrix = Math::GSL::Matrix->new(3,3);
    ...
    my $matrix_row = $matrix->row(0);

=cut

sub row
{
    my ($self, $row) = @_;
    croak (__PACKAGE__.'::$matrix->row($row) - invalid $row value')
        unless (($row < $self->rows) and $row >= 0);

   my $rowmat = Math::GSL::Matrix->new(1,$self->cols);

   for my $n (0 .. $self->cols-1) {
        gsl_matrix_set( $rowmat->raw, 0, $n,
            gsl_matrix_get($self->raw, $row, $n)
        );
   }

    return $rowmat;
}

=head3 col()

Returns a col matrix of the column you enter.

    my $matrix = Math::GSL::Matrix->new(3,3);
    ...
    my $matrix_col = $matrix->col(0);

=cut

sub col
{
    my ($self, $col) = @_;
    croak (__PACKAGE__."::\$matrix->col(\$col) - $col not a valid column")
        unless ($col < $self->cols and $col >= 0);

    my $colmat = Math::GSL::Matrix->new($self->rows, 1);

   map { gsl_matrix_set($colmat->raw, $_, 0,
           gsl_matrix_get($self->raw, $_, $col),
           )
       } (0 .. $self->rows-1);

   return $colmat;
}


=head3 max()

Computes the maximum value of a matrix. In array context returns the maximum
value and the position of that element in the matrix. If the matrix is a
vector (being it vertical or horizontal) only one coordinate is returned: the
position of the element in the vector.

    $max = $matrix->max();
    ($max, $i, $j) = $matrix->max();

=cut

sub max {
    my $matrix = shift;

    my $max = gsl_matrix_max($matrix->raw);

    if (wantarray) {
        my ($i, $j) = gsl_matrix_max_index($matrix->raw);

        return ($max, $j) if $matrix->rows == 1;
        return ($max, $i) if $matrix->cols == 1;
        return ($max, $i, $j);

    } else {
        return $max;
    }
}

=head3 min()

Computes the minimum value of a matrix. In array context returns the minimum
value and the position of that element in the matrix. If the matrix is a
vector (being it vertical or horizontal) only one coordinate is returned: the
position of the element in the vector.

    $min = $matrix->min();
    ($min, $i, $j) = $matrix->min();

=cut

sub min {
    my $matrix = shift;

    my $min = gsl_matrix_min($matrix->raw);

    if (wantarray) {
        my ($i, $j) = gsl_matrix_min_index($matrix->raw);

        return ($min, $j) if $matrix->rows == 1;
        return ($min, $i) if $matrix->cols == 1;
        return ($min, $i, $j);

    } else {
        return $min;
    }
}



=head2 Setters

=head3 identity()

Set a matrix to the identity matrix, i.e. one on the diagonal and zero elsewhere.

    my $I = $matrix->identity;

=cut

sub identity
{
    my $self=shift;
    gsl_matrix_set_identity($self->raw);
    return $self;
}

=head3 zero()

Set a matrix to the zero matrix.

    $matrix->zero;

=cut

sub zero # brrr!
{
    my $self=shift;
    gsl_matrix_set_zero($self->raw);
    return $self;
}

=head3 set_elem()

Sets a specific value in the matrix.

    my $matrix = Math::GSL::Matrix->new(2,2);
    $matrix->set_elem(0, 0, $value);

You can set multiple elements at once with chained calls:

    $matrix->set_elem(0,0,1)->set_elem(1,1,1);

B<NOTE:> just like any other method on this module, rows and arrays start
with indice 0.

=cut

sub set_elem {
    my ($self, $row, $col, $value) = @_;
    die __PACKAGE__.'::set_elem($x, $y, $v) - $x must be a valid row number'
        if ($row < 0 || $row >= $self->rows);
    die __PACKAGE__.'::set_elem($x, $y, $v) - $y must be a valid column number'
        if ($col < 0 || $col >= $self->cols);
    gsl_matrix_set($self->raw, $row, $col, $value);
    return $self;
}

=head3 set_row()

Sets a the values of a row with the elements of an array.

    my $matrix = Math::GSL::Matrix->new(3,3);
    $matrix->set_row(0, [8, 6, 2]);

You can also set multiple rows at once with chained calls:

    my $matrix = Math::GSL::Matrix->new(3,3);
    $matrix->set_row(0, [8, 6, 2])
           ->set_row(1, [2, 4, 1]);
    ...

=cut

sub set_row {
    my ($self, $row, $values) = @_;
    my $length = $#$values;
    die __PACKAGE__.'::set_row($x, $values) - $values must be a nonempty array reference' if $length == -1;
    die __PACKAGE__.'::set_row($x, $values) - $x must be a valid row number'
        if ($row < 0 || $row >= $self->rows);
    die __PACKAGE__.'::set_row($x, $values) - $values must contains the same number of elements as there is columns in the matrix'
        if($length != $self->cols-1);
    map { gsl_matrix_set($self->raw, $row, $_, $values->[$_]) } (0..$length);
    return $self;
}

=head3 set_col()

Sets a the values of a column with the elements of an array.

    my $matrix = Math::GSL::Matrix->new(3,3);
    $matrix->set_col(0, [8, 6, 2]);

You can also set multiple columns at once with chained calls:
    my $matrix = Math::GSL::Matrix->new(3,3);
    $matrix->set_col(0, [8, 6, 2])
           ->set_col(1, [2, 4, 1]);
    ...

=cut

sub set_col {
    my ($self, $col, $values) = @_;
    my $length = $#$values;

    die __PACKAGE__.'::set_col($x, $values) - $values must be a nonempty array reference' if $length == -1;
    die __PACKAGE__.'::set_col($x, $values) - $x must be a valid column number'
        if ($col < 0 || $col >= $self->cols);
    die __PACKAGE__.'::set_col($x, $values) - $values must contains the same number of elements as there is rowss in the matrix'
        if($length != $self->rows-1);

    map { gsl_matrix_set($self->raw, $_, $col, $values->[$_]) } (0..$length);
    return $self;
}

=head2 Utility Functions

=head3 copy()

Returns a copy of the matrix, which has the same size and values but resides at a different location in memory.

    my $matrix = Math::GSL::Matrix->new(5,5);
    my $copy   = $matrix->copy;

=cut

sub copy {
    my $self = shift;
    my $copy = Math::GSL::Matrix->new( $self->rows, $self->cols );
    if ( gsl_matrix_memcpy($copy->raw, $self->raw) != $GSL_SUCCESS ) {
        croak "Math::GSL - error copying memory, aborting";
    }
    return $copy;
}




=head3 is_square()

Returns true if a matrix is square, i.e. it has the same number of rows as columns, false otherwise.

=cut

sub is_square($)
{
    my $self = shift;
    return ($self->rows == $self->cols) ? 1 : 0 ;
}

=head3 det()

Returns the determinant of a matrix (computed by LU decomposition) or dies if called on a non-square matrix.

    my $det = $matrix->det();

=cut

sub det($)
{
    my $self = shift;
    croak(__PACKAGE__."- determinant only exists for square matrices") unless $self->is_square;

    my $p  = Math::GSL::Permutation->new( $self->rows );
    my $LU = $self->copy;

    my ($result, $s) = gsl_linalg_LU_decomp($LU->raw, $p->raw);
    return gsl_linalg_LU_det($LU->raw, $s );
}

=head3 write

Saves B<ONE> matrix into a file using the GSL binary format.

     $matrix->save("matrix.dat");

Note that this method always overwrite the file if it exists. In case more
than one GSL object will be written to the file,  please use the C<gsl_fopen>,
C<gsl_matrix_fwrite> and C<gsl_fclose> manually. Future versions might include
helper for this task.

B<NOTE>: in order to allow the read of one of these files without needing
to know in advance the size of the matrix, an horizontal matrix with two
elements is saved in the beginning of the file with the number of rows
and columns of the saved matrix.

=cut

sub write
{
    my ($self, $filename) = @_;

    my $fh = gsl_fopen($filename, "w");
    croak __PACKAGE__."::write - error opening file $filename for writing" unless $fh;

    # create a temporary matrix with the main matrix dimensions
    my $dim = gsl_matrix_alloc(1, 2);
    my ($rows, $cols) = $self->dim;
    gsl_matrix_set($dim, 0, 0, $rows);
    gsl_matrix_set($dim, 0, 1, $cols);

    my $err = gsl_matrix_fwrite($fh, $dim);
    croak __PACKAGE__."::write - error writing matrix" if $err;

    gsl_matrix_free($dim);

    $err = gsl_matrix_fwrite($fh, $self->raw);
    croak __PACKAGE__."::write - error writing matrix" if $err;

    gsl_fclose($fh);
}

=head3 read

Loads a matrix from a GSL binary file (saved using the C<save> method,
as it stores some extra information on the file).

    my $matrix = Math::GSL::Matrix->read("matrix.dat");

=cut

sub read
{
    my ($class, $filename) = @_;
    croak __PACKAGE__."::read - $filename does not exist" unless -f $filename;

    my $fh = gsl_fopen($filename, "r");
    croak __PACKAGE__."::read - error opening file $filename for reading" unless $fh;

    my $dim = gsl_matrix_alloc(1, 2);
    my $err = gsl_matrix_fread($fh, $dim);
    croak __PACKAGE__."::read - error reading matrix" if $err;

    my $m = gsl_matrix_alloc(gsl_matrix_get($dim, 0, 0),
                             gsl_matrix_get($dim, 0, 1));
    $err = gsl_matrix_fread($fh, $m);
    croak __PACKAGE__."::read - error reading matrix" if $err;

    gsl_matrix_free($dim);

    gsl_fclose($fh);
    return Math::GSL::Matrix->new($m);
}

=head3 lndet()

Returns the natural log of the absolute value of the determinant of a matrix (computed by LU decomposition) or dies if called on a non-square matrix.

    my $lndet = $matrix->lndet();

=cut

sub lndet($)
{
    my $self = shift;
    croak(__PACKAGE__."- log determinant only exists for square matrices") unless $self->is_square;

    my $p  = Math::GSL::Permutation->new( $self->rows );
    my $LU = $self->copy;

    gsl_linalg_LU_decomp($LU->raw, $p->raw);
    return gsl_linalg_LU_lndet($LU->raw);
}

=head3 inverse()

Returns the inverse of a matrix or dies when called on a non-square matrix.

    my $inverse = $matrix->inverse;

=cut

sub inverse($)
{
    my $self = shift;
    croak(__PACKAGE__."- inverse only exists for square matrices") unless $self->is_square;

    my $p  = Math::GSL::Permutation->new( $self->rows );
    my $LU = $self->copy;
    my $inverse = $self->copy;

    # should check return status
    gsl_linalg_LU_decomp($LU->raw, $p->raw);
    gsl_linalg_LU_invert($LU->raw, $p->raw,$inverse->raw);

    return $inverse;

}

=head3 transpose()

Returns the transpose of a matrix  or dies when called on a non-square matrix.

    my $transposed = $matrix->transpose;

=cut

sub transpose($)
{
    my $self = shift;
    croak(__PACKAGE__."- transpose only exists for square matrices") unless $self->is_square;

    my $copy = $self->copy;

    # should check return status
    gsl_matrix_transpose($copy->raw);

    return $copy;
}

=head3 eigenvalues()

Computes the a matrix eigen values.

=cut


sub eigenvalues($)
{
    my $self=shift;
    my ($r,$c) = ($self->rows,$self->cols);
    croak "Math::GSL::Matrix : \$matrix->eigenvalues - \$matrix must be square" unless ($r == $c);
    my $evec   = gsl_matrix_complex_alloc($r,$c);
    my $eigen  = gsl_eigen_nonsymmv_alloc($r);
    my $vector = gsl_vector_complex_alloc($r);

    gsl_eigen_nonsymmv($self->raw, $vector, $evec, $eigen);

    my $x = gsl_vector_complex_real($vector); # vector of real components
    my $y = gsl_vector_complex_imag($vector); # vector of imaginary components

    return map {
                _as_complex( gsl_vector_get($x->{vector}, $_) ,
                             gsl_vector_get($y->{vector}, $_) )
                } (0 .. $r-1);

}

=head3 eigenpair()

=cut


sub eigenpair($)
{
    my $self=shift;
    my ($r,$c) = ($self->rows,$self->cols);
    croak "Math::GSL::Matrix : \$matrix->eigenvalues - \$matrix must be square" unless ($r == $c);
    my $evec   = Math::GSL::MatrixComplex->new($r,$c);
    my $eigen  = gsl_eigen_nonsymmv_alloc($r);
    my $vector = Math::GSL::VectorComplex->new($r);

    gsl_eigen_nonsymmv($self->raw, $vector->raw, $evec->raw, $eigen);

    my $eigenvectors = [ map { $evec->col($_)->as_vector } (0 .. $r-1) ];

    my $x = gsl_vector_complex_real($vector->raw); # vector of real components
    my $y = gsl_vector_complex_imag($vector->raw); # vector of imaginary components
    my $eigenvalues = [
                            map {
                                _as_complex( gsl_vector_get($x->{vector}, $_) ,
                                            gsl_vector_get($y->{vector}, $_) )
                             } (0 .. $r-1)
                      ];
    return ($eigenvalues, $eigenvectors);
}

sub _as_complex {
    my ($w,$v) = @_;
    my ($x,$y) = ($w,$v);
    if( ref $w eq 'Math::GSL::Complex') {
        ($x,$y) = (gsl_real($w),gsl_imag($w));
    }
    my $z      = Math::Complex->make( $x, $y);
}


=head3 ieach()

Applied a function to each element of a matrix.

    # compute exp to each element of matrix
    $matrix->ieach( sub { exp(shift) });

This method B<changes> the original matrix. See C<each> for a
non destructive method.

=cut

sub ieach($&) {
	my ($matrix, $func) = @_;

    die __PACKAGE__.'::ieach($func) - $func must be a code reference'
    	unless ref($func) eq "CODE";


	for my $i (0 .. $matrix->rows-1) {
		for my $j (0 .. $matrix->cols-1) {
			# should be faster than use ->set_elem and ->get_elem
			gsl_matrix_set($matrix->raw, $i, $j,
				           $func->(gsl_matrix_get($matrix->raw, $i, $j)));
		}
	}
	return $matrix;
}


=head3 each()

Applied a function to each element of a matrix.

    # compute exp to each element of matrix
    $matrix->each( sub { exp(shift) });

This method returns a B<new> matrix.

=cut

sub each($&) {
	my ($matrix, $func) = @_;

    die __PACKAGE__.'::each($func) - $func must be a code reference'
    	unless ref($func) eq "CODE";

    my $result = Math::GSL::Matrix->new($matrix->rows, $matrix->cols);
	for my $i (0 .. $matrix->rows-1) {
		for my $j (0 .. $matrix->cols-1) {
			# should be faster than use ->set_elem and ->get_elem
			gsl_matrix_set($result->raw, $i, $j,
				           $func->(gsl_matrix_get($matrix->raw, $i, $j)));
		}
	}
	return $result;
}

=head3 vconcat ()

Concatenates vertically a new matrix with the object matrix, where the
object matrix is above in the final matrix.
Note that both matrices need to have the same number of columns.
The method returns a B<new> matrix.

    my $final = $m1->vconcat($m2);
    # $final is $m1 on top of $m2

=cut

sub vconcat {
    my ($self, $b) = @_;

    die __PACKAGE__.'::vconcat($m) - $m must be a Math::GSL::Matrix object'
        unless blessed($b) && $b->isa('Math::GSL::Matrix');
    die __PACKAGE__.'::vconcat($m) - $self and $m should have same number of columns'
        unless $self->cols == $b->cols;

    bless {
         _matrix => gsl_matrix_vconcat($self->raw, $b->raw),
         _rows   => $self->rows + $b->rows,
         _cols   => $self->cols
    }, 'Math::GSL::Matrix'
}

=head3 hconcat ()

Concatenates horizontally a new matrix with the object matrix, where the
object matrix is at the left in the final matrix.
Note that both matrices need to have the same number of rows.
The method returns a B<new> matrix.

    my $final = $m1->hconcat($m2);
    # $final is $m1 at the left of $m2

=cut

sub hconcat {
    my ($self, $b) = @_;

    die __PACKAGE__.'::hconcat($m) - $m must be a Math::GSL::Matrix object'
        unless blessed($b) && $b->isa('Math::GSL::Matrix');
    die __PACKAGE__.'::hconcat($m) - $self and $m should have same number of rows'
        unless $self->rows == $b->rows;

    bless {
         _matrix => gsl_matrix_hconcat($self->raw, $b->raw),
         _rows   => $self->rows,
         _cols   => $self->cols + $b->cols,
    }, 'Math::GSL::Matrix'
}

=head3 random ()

Fills a matrix with random values between 0 and 1

=cut

sub random {
	my ($self) = @_;
	
	gsl_matrix_random($self->raw);
	
	return $self;
}




sub _addition {
    my ($left, $right) = @_;

    my $lcopy = $left->copy;

    if ( blessed $right && $right->isa('Math::GSL::Matrix') && blessed $left && $left->isa('Math::GSL::Matrix') ) {
        if ( $left->cols == $right->cols && $left->rows == $right->rows ) {
            gsl_matrix_add($lcopy->raw, $right->raw);
        } else {
            croak "Math::GSL - addition of matrices must be called with two objects matrices and must have the same number of columns and rows";
        }
    } else {
        gsl_matrix_add_constant($lcopy->raw, $right);
    }
    return $lcopy;
}

sub _subtract {
    my ($left, $right) = @_;

    my $lcopy = $left->copy;

    if ( blessed $right && $right->isa('Math::GSL::Matrix') && blessed $left && $left->isa('Math::GSL::Matrix') ) {
        if ( $left->cols == $right->cols && $left->rows == $right->rows ) {
            gsl_matrix_sub($lcopy->raw, $right->raw);
        } else {
            croak "Math::GSL - subtraction of matrices must be called with two objects matrices and must have the same number of columns and rows";
        }
    } else {
        gsl_matrix_add_constant($lcopy->raw, $right*-1);
    }
    return $lcopy;
}

sub _multiplication {
    my ($left,$right) = @_;
    my $lcopy = $left->copy;

    if ( blessed $right && $right->isa('Math::GSL::Matrix') ) {
        return _dot_product($left,$right);
    } else {
        gsl_matrix_scale($lcopy->raw, $right);
    }
    return $lcopy;
}

sub _equal {
    my ($left, $right) = @_;

    return is_similar( [ $left->as_list ], [ $right->as_list ] );
}

sub _not_equal {
    my ($left, $right) = @_;
    return !is_similar( [ $left->as_list ], [ $right->as_list ] );
}

sub _dot_product {
    my ($left,$right) = @_;

    croak "Math::GSL::Matrix - incompatible matrices in multiplication"
        unless ( blessed $right && $right->isa('Math::GSL::Matrix') and $left->cols == $right->rows );
    my $C = Math::GSL::Matrix->new($left->rows, $right->cols);
    gsl_blas_dgemm($CblasNoTrans, $CblasNoTrans, 1, $left->raw, $right->raw, 0, $C->raw);
    return $C;
}

sub DESTROY {
    my $self = shift;

    # during global phase of destruction there is no
    # sure that objects are destroyed in order. So, when
    # going out, just do not destroy anything. The Operating
    # System should do us for us.
    return if defined(${^GLOBAL_PHASE}) && ${^GLOBAL_PHASE} eq 'DESTRUCT';

    gsl_matrix_free($self->raw) if defined($self->raw);
}

=head1 GSL FUNCTION INTERFACE

Here is a list of all the functions included in this module :

=over 1

=item C<gsl_matrix_alloc($i, $j)> - Return a gsl_matrix of $i rows and $j columns

=item C<gsl_matrix_calloc($i, $j)> - Return a gsl_matrix of $i rows and $j columns and initialize all of the elements of the matrix to zero

=item C<gsl_matrix_alloc_from_block> -

=item C<gsl_matrix_free> -

=item C<gsl_matrix_alloc_from_matrix > -

=item C<gsl_vector_alloc_row_from_matrix> -

=item C<gsl_vector_alloc_col_from_matrix > -

=item C<gsl_matrix_submatrix($m, $k1, $k2, $n1, $n2)> - Return a matrix view of the matrix $m. The upper-left element of the submatrix is the element ($k1,$k2) of the original matrix. The submatrix has $n1 rows and $n2 columns.

=item C<gsl_matrix_row($m , $i)> - Return a vector view of the $i-th row of the matrix $m

=item C<gsl_matrix_column($m, $j)> - Return a vector view of the $j-th column of the matrix $m

=item C<gsl_matrix_diagonal($m)> - Return a vector view of the diagonal of the vector. The matrix doesn't have to be square.

=item C<gsl_matrix_subdiagonal($m, $k)> - Return a vector view of the $k-th subdiagonal of the matrix $m. The diagonal of the matrix corresponds to k=0.

=item C<gsl_matrix_superdiagonal($m, $k)> - Return a vector view of the $k-th superdiagonal of the matrix $m. The matrix doesn't have to be square.

=item C<gsl_matrix_subrow($m, $i, $offset, $n)> - Return a vector view of the $i-th row of the matrix $m beginning at offset elements and containing n elements.

=item C<gsl_matrix_subcolumn($m, $j, $offset, $n)> - Return a vector view of the $j-th column of the matrix $m beginning at offset elements and containing n elements.

=item C<gsl_matrix_view_array($base, $n1, $n2)> - This function returns a matrix view of the array reference $base. The matrix has $n1 rows and $n2 columns. The physical number of columns in memory is also given by $n2. Mathematically, the (i,j)-th element of the new matrix is given by, m'(i,j) = $base->[i*$n2 + j] where the index i runs from 0 to $n1-1 and the index j runs from 0 to $n2-1. The new matrix is only a view of the array reference $base. When the view goes out of scope the original array reference $base will continue to exist. The original memory can only be deallocated by freeing the original array. Of course, the original array should not be deallocated while the view is still in use.

=item C<gsl_matrix_view_array_with_tda($base, $n1, $n2, $tda)> - This function returns a matrix view of the array reference $base with a physical number of columns $tda which may differ from the corresponding dimension of the matrix. The matrix has $n1 rows and $n2 columns, and the physical number of columns in memory is given by $tda. Mathematically, the (i,j)-th element of the new matrix is given by, m'(i,j) = $base->[i*$tda + j] where the index i runs from 0 to $n1-1 and the index j runs from 0 to $n2-1. The new matrix is only a view of the array reference $base. When the view goes out of scope the original array reference $base will continue to exist. The original memory can only be deallocated by freeing the original array. Of course, the original array should not be deallocated while the view is still in use.

=item C<gsl_matrix_view_vector> -

=item C<gsl_matrix_view_vector_with_tda> -

=item C<gsl_matrix_const_submatrix> -

=item C<gsl_matrix_get($m, $i, $j)> - Return the (i,j)-th element of the matrix $m

=item C<gsl_matrix_set($m, $i, $j, $x)> - Set the value of the (i,j)-th element of the matrix $m to $x

=item C<gsl_matrix_ptr> -

=item C<gsl_matrix_const_ptr> -

=item C<gsl_matrix_set_zero($m)> - Set all the elements of the matrix $m to zero

=item C<gsl_matrix_set_identity($m)> - Set the elements of the matrix $m to the corresponding elements of the identity matrix

=item C<gsl_matrix_set_all($m, $x)> - Set all the elements of the matrix $m to the value $x

=item C<gsl_matrix_fread($fh, $m)> - Read a file which has been written with gsl_matrix_fwrite from the stream $fh opened with the gsl_fopen function from the Math::GSL module and stores the data inside the matrix $m

=item C<gsl_matrix_fwrite($fh, $m)> - Write the elements of the matrix $m in binary format to a stream $fh opened with the gsl_fopen function from the Math::GSL module

=item C<gsl_matrix_fscanf($fh, $m)> - Read a file which has been written with gsl_matrix_fprintf from the stream $fh opened with the gsl_fopenfunction from the Math::GSL module and stores the data inside the matrix $m

=item C<gsl_matrix_fprintf($fh, $m, $format)> - Write the elements of the matrix $m in the format $format (for example "%f" is the format for double) to a stream $fh opened with the gsl_fopen function from the Math::GSL module

=item C<gsl_matrix_memcpy($dest, $src)> - Copy the elements of the matrix $src to the matrix $dest. The two matrices must have the same size.

=item C<gsl_matrix_swap($m1, $m2)> - Exchange the elements of the matrices $m1 and $m2 by copying. The two matrices must have the same size.

=item C<gsl_matrix_swap_rows($m, $i, $j)> - Exchange the $i-th and $j-th row of the matrix $m. The function returns 0 if the operation succeeded, 1 otherwise.

=item C<gsl_matrix_swap_columns($m, $i, $j)> - Exchange the $i-th and $j-th column of the matrix $m. The function returns 0 if the operation succeeded, 1 otherwise.

=item C<gsl_matrix_swap_rowcol($m, $i, $j)> - Exchange the $i-th row and the $j-th column of the matrix $m. The matrix must be square. The function returns 0 if the operation succeeded, 1 otherwise.

=item C<gsl_matrix_transpose($m)> - This function replaces the matrix m by its transpose by copying the elements of the matrix in-place. The matrix must be square for this operation to be possible.

=item C<gsl_matrix_transpose_memcpy($dest, $src)> - Make the matrix $dest the transpose of the matrix $src. This function works for all matrices provided that the dimensions of the matrix dest match the transposed dimensions of the matrix src.

=item C<gsl_matrix_max($m)> - Return the maximum value in the matrix $m

=item C<gsl_matrix_min($m)> - Return the minimum value in the matrix $m

=item C<gsl_matrix_minmax($m)> - Return two values, the first is the minimum value of the Matrix $m and the second is the maximum of the same the same matrix.

=item C<gsl_matrix_max_index($m)> - Return two values, the first is the the i indice of the maximum value of the matrix $m and the second is the j indice of the same value.

=item C<gsl_matrix_min_index($m)> - Return two values, the first is the the i indice of the minimum value of the matrix $m and the second is the j indice of the same value.

=item C<gsl_matrix_minmax_index($m)> - Return four values, the first is the i indice of the minimum of the matrix $m, the second is the j indice of the same value, the third is the i indice of the maximum of the matrix $m and the fourth is the j indice of the same value

=item C<gsl_matrix_isnull($m)> - Return 1 if all the elements of the matrix $m are zero, 0 otherwise

=item C<gsl_matrix_ispos($m)> - Return 1 if all the elements of the matrix $m are strictly positive, 0 otherwise

=item C<gsl_matrix_isneg($m)> - Return 1 if all the elements of the matrix $m are strictly negative, 0 otherwise

=item C<gsl_matrix_isnonneg($m)> - Return 1 if all the elements of the matrix $m are non-negatuive, 0 otherwise

=item C<gsl_matrix_add($a, $b)> - Add the elements of matrix $b to the elements of matrix $a

=item C<gsl_matrix_sub($a, $b)> - Subtract the elements of matrix $b from the elements of matrix $a

=item C<gsl_matrix_mul_elements($a, $b)> - Multiplie the elements of matrix $a by the elements of matrix $b

=item C<gsl_matrix_div_elements($a, $b)> - Divide the elements of matrix $a by the elements of matrix $b

=item C<gsl_matrix_scale($a, $x)> - Multiplie the elements of matrix $a by the constant factor $x

=item C<gsl_matrix_add_constant($a, $x)> - Add the constant value $x to the elements of the matrix $a

=item C<gsl_matrix_add_diagonal($a, $x)> - Add the constant value $x to the elements of the diagonal of the matrix $a

=item C<gsl_matrix_get_row($v, $m, $i)> - Copy the elements of the $i-th row of the matrix $m into the vector $v. The length of the vector must be of the same as the length of the row. The function returns 0 if it succeeded, 1 otherwise.

=item C<gsl_matrix_get_col($v, $m, $i)> - Copy the elements of the $j-th column of the matrix $m into the vector $v. The length of the vector must be of the same as the length of the column. The function returns 0 if it succeeded, 1 otherwise.

=item C<gsl_matrix_set_row($m, $i, $v)> - Copy the elements of vector $v into the $i-th row of the matrix $m The length of the vector must be of the same as the length of the row. The function returns 0 if it succeeded, 1 otherwise.

=item C<gsl_matrix_set_col($m, $j, $v)> - Copy the elements of vector $v into the $j-th row of the matrix $m The length of the vector must be of the same as the length of the column. The function returns 0 if it succeeded, 1 otherwise.

=back

These are related to constant views of a matrix.

=over 1


=item C<gsl_matrix_const_row>

=item C<gsl_matrix_const_column>

=item C<gsl_matrix_const_diagonal>

=item C<gsl_matrix_const_subdiagonal>

=item C<gsl_matrix_const_superdiagonal>

=item C<gsl_matrix_const_subrow>

=item C<gsl_matrix_const_subcolumn>

=item C<gsl_matrix_const_view_array>

=item C<gsl_matrix_const_view_array_with_tda>

=back


The following functions are similar to those above but work with C<char>'s and C<int>'s. We are not quite
sure if anyone wants these. Please speak up if you do and/or submit some patches to this documentation, please!

=over 1


=item gsl_matrix_const_view_vector

=item gsl_matrix_const_view_vector_with_tda

=item gsl_matrix_char_alloc

=item gsl_matrix_char_calloc

=item gsl_matrix_char_alloc_from_block

=item gsl_matrix_char_alloc_from_matrix

=item gsl_vector_char_alloc_row_from_matrix

=item gsl_vector_char_alloc_col_from_matrix

=item gsl_matrix_char_free

=item gsl_matrix_char_submatrix

=item gsl_matrix_char_row

=item gsl_matrix_char_column

=item gsl_matrix_char_diagonal

=item gsl_matrix_char_subdiagonal

=item gsl_matrix_char_superdiagonal

=item gsl_matrix_char_subrow

=item gsl_matrix_char_subcolumn

=item gsl_matrix_char_view_array

=item gsl_matrix_char_view_array_with_tda

=item gsl_matrix_char_view_vector

=item gsl_matrix_char_view_vector_with_tda

=item gsl_matrix_char_const_submatrix

=item gsl_matrix_char_const_row

=item gsl_matrix_char_const_column

=item gsl_matrix_char_const_diagonal

=item gsl_matrix_char_const_subdiagonal

=item gsl_matrix_char_const_superdiagonal

=item gsl_matrix_char_const_subrow

=item gsl_matrix_char_const_subcolumn

=item gsl_matrix_char_const_view_array

=item gsl_matrix_char_const_view_array_with_tda

=item gsl_matrix_char_const_view_vector

=item gsl_matrix_char_const_view_vector_with_tda

=item gsl_matrix_char_get

=item gsl_matrix_char_set

=item gsl_matrix_char_ptr

=item gsl_matrix_char_const_ptr

=item gsl_matrix_char_set_zero

=item gsl_matrix_char_set_identity

=item gsl_matrix_char_set_all

=item gsl_matrix_char_fread

=item gsl_matrix_char_fwrite

=item gsl_matrix_char_fscanf

=item gsl_matrix_char_fprintf

=item gsl_matrix_char_memcpy

=item gsl_matrix_char_swap

=item gsl_matrix_char_swap_rows

=item gsl_matrix_char_swap_columns

=item gsl_matrix_char_swap_rowcol

=item gsl_matrix_char_transpose

=item gsl_matrix_char_transpose_memcpy

=item gsl_matrix_char_max

=item gsl_matrix_char_min

=item gsl_matrix_char_minmax

=item gsl_matrix_char_max_index

=item gsl_matrix_char_min_index

=item gsl_matrix_char_minmax_index

=item gsl_matrix_char_isnull

=item gsl_matrix_char_ispos

=item gsl_matrix_char_isneg

=item gsl_matrix_char_isnonneg

=item gsl_matrix_char_add

=item gsl_matrix_char_sub

=item gsl_matrix_char_mul_elements

=item gsl_matrix_char_div_elements

=item gsl_matrix_char_scale

=item gsl_matrix_char_add_constant

=item gsl_matrix_char_add_diagonal

=item gsl_matrix_int_alloc

=item gsl_matrix_int_calloc

=item gsl_matrix_int_alloc_from_block

=item gsl_matrix_int_alloc_from_matrix

=item gsl_vector_int_alloc_row_from_matrix

=item gsl_vector_int_alloc_col_from_matrix

=item gsl_matrix_int_free

=item gsl_matrix_int_submatrix

=item gsl_matrix_int_row

=item gsl_matrix_int_column

=item gsl_matrix_int_diagonal

=item gsl_matrix_int_subdiagonal

=item gsl_matrix_int_superdiagonal

=item gsl_matrix_int_subrow

=item gsl_matrix_int_subcolumn

=item gsl_matrix_int_view_array

=item gsl_matrix_int_view_array_with_tda

=item gsl_matrix_int_view_vector

=item gsl_matrix_int_view_vector_with_tda

=item gsl_matrix_int_const_submatrix

=item gsl_matrix_int_const_row

=item gsl_matrix_int_const_column

=item gsl_matrix_int_ptr

=item gsl_matrix_int_const_ptr

=item gsl_matrix_int_set_zero

=item gsl_matrix_int_set_identity

=item gsl_matrix_int_set_all

=item gsl_matrix_int_fread

=item gsl_matrix_int_fwrite

=item gsl_matrix_int_fscanf

=item gsl_matrix_int_fprintf

=item gsl_matrix_int_memcpy

=item gsl_matrix_int_swap

=item gsl_matrix_int_swap_rows

=item gsl_matrix_int_swap_columns

=item gsl_matrix_int_swap_rowcol

=item gsl_matrix_int_transpose

=item gsl_matrix_int_transpose_memcpy

=item gsl_matrix_int_max

=item gsl_matrix_int_min

=item gsl_matrix_int_minmax

=item gsl_matrix_int_max_index

=item gsl_matrix_int_min_index

=item gsl_matrix_int_minmax_index

=item gsl_matrix_int_isnull

=item gsl_matrix_int_ispos

=item gsl_matrix_int_isneg

=item gsl_matrix_int_isnonneg

=item gsl_matrix_int_add

=item gsl_matrix_int_sub

=item gsl_matrix_int_mul_elements

=item gsl_matrix_int_div_elements

=item gsl_matrix_int_scale

=item gsl_matrix_int_add_constant

=item gsl_matrix_int_add_diagonal

=back

You have to add the functions you want to use inside the qw /put_function_here /.
You can also write use Math::GSL::Matrix qw/:all/ to use all available functions of the module.
Other tags are also available, here is a complete list of all tags for this module :

=over 1

=item C<all>

=item C<int>

=item C<double>

=item C<char>

=item C<complex>

=back

For more information on the functions, we refer you to the GSL official documentation
L<http://www.gnu.org/software/gsl/manual/html_node/>



=head1 EXAMPLES

 Most of the examples from this section are perl versions of the examples at L<http://www.gnu.org/software/gsl/manual/html_node/Example-programs-for-matrices.html>

 The program below shows how to allocate, initialize and read from a matrix using the functions gsl_matrix_alloc, gsl_matrix_set and gsl_matrix_get.

 use Math::GSL::Matrix qw/:all/;
 my $m = gsl_matrix_alloc (10,3);
 for my $i (0..9){
    for my $j (0..2){
        gsl_matrix_set($m, $i, $j, 0.23 + 100*$i + $j);
    }
 }

 for my $i (0..99){ # OUT OF RANGE ERROR
     for my $j (0..2){
        print "m($i, $j) = " . gsl_matrix_get ($m, $i, $j) . "\n";
    }
 }
 gsl_matrix_free ($m);


 use Math::GSL::Matrix qw/:all/;

 my $m = gsl_matrix_alloc (100, 100);
 my $a = gsl_matrix_alloc (100, 100);

 for my $i (0..99){
     for my $j (0..99){
         gsl_matrix_set ($m, $i, $j, 0.23 + $i + $j);
     }
 }

 The next program shows how to write a matrix to a file.

 my $out = gsl_fopen("test.dat", "wb");
 gsl_matrix_fwrite ($out, $m);
 gsl_fclose ($out);

 my $in = gsl_fopen("test.dat", "rb");
 gsl_matrix_fread ($in, $a);
 gsl_fclose($in);

 my $k=0;
 for my $i (0..99){
     for my $j (0..99){
         $mij = gsl_matrix_get ($m, $i, $j);
         $aij = gsl_matrix_get ($a, $i, $j);
         $k++ if ($mij != $aij);
     }
 }

 gsl_matrix_free($m);
 gsl_matrix_free($a);

 print "differences = $k (should be zero)\n";

=head1 AUTHORS

Jonathan "Duke" Leto <jonathan@leto.net> and Thierry Moisan <thierry.moisan@gmail.com>

=head1 COPYRIGHT AND LICENSE

Copyright (C) 2008-2024 Jonathan "Duke" Leto and Thierry Moisan

This program is free software; you can redistribute it and/or modify it
under the same terms as Perl itself.

=cut

%}