6 use Scalar::Util 'blessed';
8 use Math::GSL qw/:all/;
9 use Math::GSL::Errno qw/:all/;
10 use Math::GSL::Eigen qw/:all/;
11 use Math::GSL::Vector qw/:all/;
12 use Math::GSL::Complex qw/:all/;
13 use Math::GSL::Test qw/is_similar/;
14 use Math::GSL::BLAS qw/gsl_blas_dgemm/;
15 use Math::GSL::CBLAS qw/$CblasNoTrans/;
16 use Math::GSL::Permutation;
17 use Math::GSL::Linalg qw/
24 # should only include needed methods
25 use Math::GSL::MatrixComplex qw/:all/;
26 use Math::GSL::VectorComplex qw/:all/;
30 '*' => \&_multiplication,
38 gsl_matrix_alloc gsl_matrix_calloc gsl_matrix_alloc_from_block
39 gsl_matrix_alloc_from_matrix gsl_vector_alloc_row_from_matrix
40 gsl_vector_alloc_col_from_matrix gsl_matrix_free gsl_matrix_submatrix
41 gsl_matrix_row gsl_matrix_column gsl_matrix_diagonal
42 gsl_matrix_subdiagonal gsl_matrix_superdiagonal gsl_matrix_subrow
43 gsl_matrix_subcolumn gsl_matrix_view_array
44 gsl_matrix_view_array_with_tda gsl_matrix_view_vector
45 gsl_matrix_view_vector_with_tda gsl_matrix_const_submatrix
46 gsl_matrix_const_row gsl_matrix_const_column gsl_matrix_const_diagonal
47 gsl_matrix_const_subdiagonal gsl_matrix_const_superdiagonal
48 gsl_matrix_const_subrow gsl_matrix_const_subcolumn
49 gsl_matrix_const_view_array gsl_matrix_const_view_array_with_tda
50 gsl_matrix_const_view_vector gsl_matrix_const_view_vector_with_tda
51 gsl_matrix_get gsl_matrix_set gsl_matrix_ptr gsl_matrix_const_ptr
52 gsl_matrix_set_zero gsl_matrix_set_identity gsl_matrix_set_all
53 gsl_matrix_fread gsl_matrix_fwrite gsl_matrix_fscanf gsl_matrix_fprintf
54 gsl_matrix_memcpy gsl_matrix_swap gsl_matrix_swap_rows
55 gsl_matrix_swap_columns gsl_matrix_swap_rowcol gsl_matrix_transpose
56 gsl_matrix_transpose_memcpy gsl_matrix_max gsl_matrix_min
57 gsl_matrix_minmax gsl_matrix_max_index gsl_matrix_min_index
58 gsl_matrix_minmax_index gsl_matrix_isnull gsl_matrix_ispos
59 gsl_matrix_isneg gsl_matrix_isnonneg gsl_matrix_add gsl_matrix_sub
60 gsl_matrix_mul_elements gsl_matrix_div_elements gsl_matrix_scale
61 gsl_matrix_add_constant gsl_matrix_add_diagonal
62 gsl_matrix_char_alloc gsl_matrix_char_calloc gsl_matrix_char_alloc_from_block
63 gsl_matrix_char_alloc_from_matrix gsl_vector_char_alloc_row_from_matrix gsl_vector_char_alloc_col_from_matrix
64 gsl_matrix_char_free gsl_matrix_char_submatrix
65 gsl_matrix_char_row gsl_matrix_char_column
66 gsl_matrix_char_diagonal gsl_matrix_char_subdiagonal gsl_matrix_char_superdiagonal
67 gsl_matrix_char_subrow gsl_matrix_char_subcolumn gsl_matrix_char_view_array
68 gsl_matrix_char_view_array_with_tda gsl_matrix_char_view_vector gsl_matrix_char_view_vector_with_tda
69 gsl_matrix_char_const_submatrix gsl_matrix_char_const_row gsl_matrix_char_const_column
70 gsl_matrix_char_const_diagonal gsl_matrix_char_const_subdiagonal gsl_matrix_char_const_superdiagonal
71 gsl_matrix_char_const_subrow gsl_matrix_char_const_subcolumn gsl_matrix_char_const_view_array
72 gsl_matrix_char_const_view_array_with_tda gsl_matrix_char_const_view_vector gsl_matrix_char_const_view_vector_with_tda
73 gsl_matrix_char_get gsl_matrix_char_set gsl_matrix_char_ptr gsl_matrix_char_const_ptr
74 gsl_matrix_char_set_zero gsl_matrix_char_set_identity
75 gsl_matrix_char_set_all gsl_matrix_char_fread
76 gsl_matrix_char_fwrite gsl_matrix_char_fscanf gsl_matrix_char_fprintf
77 gsl_matrix_char_memcpy gsl_matrix_char_swap
78 gsl_matrix_char_swap_rows gsl_matrix_char_swap_columns
79 gsl_matrix_char_swap_rowcol gsl_matrix_char_transpose gsl_matrix_char_transpose_memcpy
80 gsl_matrix_char_max gsl_matrix_char_min
81 gsl_matrix_char_minmax gsl_matrix_char_max_index
82 gsl_matrix_char_min_index gsl_matrix_char_minmax_index
83 gsl_matrix_char_isnull gsl_matrix_char_ispos gsl_matrix_char_isneg
84 gsl_matrix_char_isnonneg gsl_matrix_char_add
85 gsl_matrix_char_sub gsl_matrix_char_mul_elements gsl_matrix_char_div_elements
86 gsl_matrix_char_scale gsl_matrix_char_add_constant gsl_matrix_char_add_diagonal
87 gsl_matrix_int_alloc gsl_matrix_int_calloc gsl_matrix_int_alloc_from_block
88 gsl_matrix_int_alloc_from_matrix gsl_vector_int_alloc_row_from_matrix gsl_vector_int_alloc_col_from_matrix
89 gsl_matrix_int_free gsl_matrix_int_submatrix gsl_matrix_int_row
90 gsl_matrix_int_column gsl_matrix_int_diagonal gsl_matrix_int_subdiagonal
91 gsl_matrix_int_superdiagonal gsl_matrix_int_subrow gsl_matrix_int_subcolumn gsl_matrix_int_view_array
92 gsl_matrix_int_view_array_with_tda gsl_matrix_int_view_vector gsl_matrix_int_view_vector_with_tda
93 gsl_matrix_int_const_submatrix gsl_matrix_int_const_row gsl_matrix_int_const_column
94 gsl_matrix_int_const_diagonal gsl_matrix_int_const_subdiagonal gsl_matrix_int_const_superdiagonal
95 gsl_matrix_int_const_subrow gsl_matrix_int_const_subcolumn gsl_matrix_int_const_view_array
96 gsl_matrix_int_const_view_array_with_tda gsl_matrix_int_const_view_vector gsl_matrix_int_const_view_vector_with_tda
97 gsl_matrix_int_get gsl_matrix_int_set
98 gsl_matrix_int_ptr gsl_matrix_int_const_ptr
99 gsl_matrix_int_set_zero gsl_matrix_int_set_identity gsl_matrix_int_set_all
100 gsl_matrix_int_fread gsl_matrix_int_fwrite
101 gsl_matrix_int_fscanf gsl_matrix_int_fprintf
102 gsl_matrix_int_memcpy gsl_matrix_int_swap
103 gsl_matrix_int_swap_rows gsl_matrix_int_swap_columns gsl_matrix_int_swap_rowcol
104 gsl_matrix_int_transpose gsl_matrix_int_transpose_memcpy
105 gsl_matrix_int_max gsl_matrix_int_min gsl_matrix_int_minmax
106 gsl_matrix_int_max_index gsl_matrix_int_min_index
107 gsl_matrix_int_minmax_index gsl_matrix_int_isnull
108 gsl_matrix_int_ispos gsl_matrix_int_isneg gsl_matrix_int_isnonneg
109 gsl_matrix_int_add gsl_matrix_int_sub
110 gsl_matrix_int_mul_elements gsl_matrix_int_div_elements gsl_matrix_int_scale
111 gsl_matrix_int_add_constant gsl_matrix_int_add_diagonal
112 gsl_matrix_get_row gsl_matrix_get_col gsl_matrix_set_row gsl_matrix_set_col
115 our %EXPORT_TAGS = ( all => [ @EXPORT_OK ],
117 gsl_matrix_char_alloc
118 gsl_matrix_char_calloc
119 gsl_matrix_char_alloc_from_block
120 gsl_matrix_char_alloc_from_matrix
121 gsl_vector_char_alloc_row_from_matrix
122 gsl_vector_char_alloc_col_from_matrix
124 gsl_matrix_char_submatrix
126 gsl_matrix_char_column
127 gsl_matrix_char_diagonal
128 gsl_matrix_char_subdiagonal
129 gsl_matrix_char_superdiagonal
130 gsl_matrix_char_subrow
131 gsl_matrix_char_subcolumn
132 gsl_matrix_char_view_array
133 gsl_matrix_char_view_array_with_tda
134 gsl_matrix_char_view_vector
135 gsl_matrix_char_view_vector_with_tda
136 gsl_matrix_char_const_submatrix
137 gsl_matrix_char_const_row
138 gsl_matrix_char_const_column
139 gsl_matrix_char_const_diagonal
140 gsl_matrix_char_const_subdiagonal
141 gsl_matrix_char_const_superdiagonal
142 gsl_matrix_char_const_subrow
143 gsl_matrix_char_const_subcolumn
144 gsl_matrix_char_const_view_array
145 gsl_matrix_char_const_view_array_with_tda
146 gsl_matrix_char_const_view_vector
147 gsl_matrix_char_const_view_vector_with_tda
151 gsl_matrix_char_const_ptr
152 gsl_matrix_char_set_zero
153 gsl_matrix_char_set_identity
154 gsl_matrix_char_set_all
155 gsl_matrix_char_fread
156 gsl_matrix_char_fwrite
157 gsl_matrix_char_fscanf
158 gsl_matrix_char_fprintf
159 gsl_matrix_char_memcpy
161 gsl_matrix_char_swap_rows
162 gsl_matrix_char_swap_columns
163 gsl_matrix_char_swap_rowcol
164 gsl_matrix_char_transpose
165 gsl_matrix_char_transpose_memcpy
168 gsl_matrix_char_minmax
169 gsl_matrix_char_max_index
170 gsl_matrix_char_min_index
171 gsl_matrix_char_minmax_index
172 gsl_matrix_char_isnull
173 gsl_matrix_char_ispos
174 gsl_matrix_char_isneg
175 gsl_matrix_char_isnonneg
178 gsl_matrix_char_mul_elements
179 gsl_matrix_char_div_elements
180 gsl_matrix_char_scale
181 gsl_matrix_char_add_constant
182 gsl_matrix_char_add_diagonal
188 gsl_matrix_alloc_from_block
189 gsl_matrix_alloc_from_matrix
190 gsl_vector_alloc_row_from_matrix
191 gsl_vector_alloc_col_from_matrix
197 gsl_matrix_subdiagonal
198 gsl_matrix_superdiagonal
201 gsl_matrix_view_array
202 gsl_matrix_view_array_with_tda
203 gsl_matrix_view_vector
204 gsl_matrix_view_vector_with_tda
205 gsl_matrix_const_submatrix
207 gsl_matrix_const_column
208 gsl_matrix_const_diagonal
209 gsl_matrix_const_subdiagonal
210 gsl_matrix_const_superdiagonal
211 gsl_matrix_const_subrow
212 gsl_matrix_const_subcolumn
213 gsl_matrix_const_view_array
214 gsl_matrix_const_view_array_with_tda
215 gsl_matrix_const_view_vector
216 gsl_matrix_const_view_vector_with_tda
222 gsl_matrix_set_identity
231 gsl_matrix_swap_columns
232 gsl_matrix_swap_rowcol
234 gsl_matrix_transpose_memcpy
239 gsl_matrix_minmax_index
245 gsl_matrix_mul_elements
246 gsl_matrix_div_elements
248 gsl_matrix_add_constant
249 gsl_matrix_add_diagonal
253 gsl_matrix_int_alloc_from_matrix
255 gsl_matrix_int_column
256 gsl_matrix_int_superdiagonal
257 gsl_matrix_int_view_array_with_tda
258 gsl_matrix_int_const_submatrix
259 gsl_matrix_int_const_diagonal
260 gsl_matrix_int_const_subrow
261 gsl_matrix_int_const_view_array_with_tda
264 gsl_matrix_int_set_zero
266 gsl_matrix_int_fscanf
267 gsl_matrix_int_memcpy
268 gsl_matrix_int_swap_rows
269 gsl_matrix_int_transpose
271 gsl_matrix_int_max_index
272 gsl_matrix_int_minmax_index
275 gsl_matrix_int_mul_elements
276 gsl_matrix_int_add_constant
283 Math::GSL::Matrix - Mathematical functions concerning Matrices
287 use Math::GSL::Matrix qw/:all/;
288 my $matrix1 = Math::GSL::Matrix->new(5,5); # OO interface
289 my $matrix2 = $matrix1 + 4; # You can add or substract values or matrices to OO matrices
290 my $matrix3 = $matrix1 - 4;
291 my $matrix4 = $matrix2 + $matrix1;
292 my $matrix5 = $matrix2 . $matrix1; # This is a scalar product, it simply multiply each element
293 # with the element of $matrix1 that have the same position
294 # See Math::GSL::BLAS if you want scalar product
296 my $matrix6 = $matrix2 . 8; # Multiply every elements of $matrix2 by 8
297 my $matrix7 = $matrix2 * $matrix1; # scalar product of two matrices
298 if($matrix1 == $matrix4) ...
299 if($matrix1 != $matrix3) ...
300 my $matrix8 = gsl_matrix_alloc(5,5); # standard interface
303 =head1 Objected Oriented Interface to GSL Math::GSL::Matrix
305 =head2 Math::GSL::Matrix->new()
307 Creates a new Matrix of the given size.
309 my $matrix = Math::GSL::Matrix->new(10,10);
314 my ($class, $rows, $cols) = @_;
317 if ( defined $rows && defined $cols &&
318 $rows > 0 && $cols > 0 &&
319 (int $rows == $rows) && (int $cols == $cols)){
321 $matrix = gsl_matrix_alloc($rows,$cols);
323 croak( __PACKAGE__.'::new($x,$y) - $x and $y must be positive integers');
325 gsl_matrix_set_zero($matrix);
326 $this->{_matrix} = $matrix;
327 ($this->{_rows}, $this->{_cols}) = ($rows,$cols);
332 Get the underlying GSL matrix object created by SWIG, useful for using gsl_matrix_* functions which do not have an OO counterpart.
334 my $matrix = Math::GSL::Matrix->new(3,3);
335 my $gsl_matrix = $matrix->raw;
336 my $stuff = gsl_matrix_get($gsl_matrix, 1, 2);
339 sub raw { (shift)->{_matrix} }
343 Set a matrix to the identity matrix, i.e. one on the diagonal and zero elsewhere.
345 my $I = $matrix->identity;
352 gsl_matrix_set_identity($self->raw);
358 Set a matrix to the zero matrix.
367 gsl_matrix_set_zero($self->raw);
373 Returns a copy of the matrix, which has the same size and values but resides at a different location in memory.
375 my $matrix = Math::GSL::Matrix->new(5,5);
376 my $copy = $matrix->copy;
382 my $copy = Math::GSL::Matrix->new( $self->rows, $self->cols );
383 if ( gsl_matrix_memcpy($copy->raw, $self->raw) != $GSL_SUCCESS ) {
384 croak "Math::GSL - error copying memory, aborting";
392 Returns the number of rows in the matrix.
394 my $rows = $matrix->rows;
397 sub rows { (shift)->{_rows} }
401 Returns the number of columns in the matrix.
403 my $cols = $matrix->cols;
406 sub cols { (shift)->{_cols} }
410 Get the contents of a Math::GSL::Matrix object as a Perl list.
412 my $matrix = Math::GSL::Matrix->new(3,3);
414 my @matrix = $matrix->as_list;
420 Returns true if a matrix is square, i.e. it has the same number of rows as columns, false otherwise.
427 return ($self->rows == $self->cols) ? 1 : 0 ;
432 Returns the determinant of a matrix (computed by LU decomposition) or dies if called on a non-square matrix.
434 my $det = $matrix->det();
441 croak(__PACKAGE__."- determinant only exists for square matrices") unless $self->is_square;
443 my $p = Math::GSL::Permutation->new( $self->rows );
444 my $LU = $self->copy;
446 my ($result, $s) = gsl_linalg_LU_decomp($LU->raw, $p->raw);
447 return gsl_linalg_LU_det($LU->raw, $s );
452 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.
454 my $lndet = $matrix->lndet();
461 croak(__PACKAGE__."- log determinant only exists for square matrices") unless $self->is_square;
463 my $p = Math::GSL::Permutation->new( $self->rows );
464 my $LU = $self->copy;
466 gsl_linalg_LU_decomp($LU->raw, $p->raw);
467 return gsl_linalg_LU_lndet($LU->raw);
472 Returns the inverse of a matrix or dies when called on a non-square matrix.
474 my $inverse = $matrix->inverse;
481 croak(__PACKAGE__."- inverse only exists for square matrices") unless $self->is_square;
483 my $p = Math::GSL::Permutation->new( $self->rows );
484 my $LU = $self->copy;
485 my $inverse = $self->copy;
487 # should check return status
488 gsl_linalg_LU_decomp($LU->raw, $p->raw);
489 gsl_linalg_LU_invert($LU->raw, $p->raw,$inverse->raw);
503 my ($r,$c) = ($self->rows,$self->cols);
504 croak "Math::GSL::Matrix : \$matrix->eigenvalues - \$matrix must be square" unless ($r == $c);
505 my $evec = gsl_matrix_complex_alloc($r,$c);
506 my $eigen = gsl_eigen_nonsymmv_alloc($r);
507 my $vector = gsl_vector_complex_alloc($r);
509 gsl_eigen_nonsymmv($self->raw, $vector, $evec, $eigen);
511 my $x = gsl_vector_complex_real($vector); # vector of real components
512 my $y = gsl_vector_complex_imag($vector); # vector of imaginary components
515 _as_complex( gsl_vector_get($x->{vector}, $_) ,
516 gsl_vector_get($y->{vector}, $_) )
524 my ($r,$c) = ($self->rows,$self->cols);
525 croak "Math::GSL::Matrix : \$matrix->eigenvalues - \$matrix must be square" unless ($r == $c);
526 my $evec = Math::GSL::MatrixComplex->new($r,$c);
527 my $eigen = gsl_eigen_nonsymmv_alloc($r);
528 my $vector = Math::GSL::VectorComplex->new($r);
530 gsl_eigen_nonsymmv($self->raw, $vector->raw, $evec->raw, $eigen);
532 my $eigenvectors = [ map { $evec->col($_)->as_vector } (0 .. $r-1) ];
534 my $x = gsl_vector_complex_real($vector->raw); # vector of real components
535 my $y = gsl_vector_complex_imag($vector->raw); # vector of imaginary components
538 _as_complex( gsl_vector_get($x->{vector}, $_) ,
539 gsl_vector_get($y->{vector}, $_) )
542 return ($eigenvalues, $eigenvectors);
547 my ($x,$y) = ($w,$v);
548 if( ref $w eq 'Math::GSL::Complex') {
549 ($x,$y) = (gsl_real($w),gsl_imag($w));
551 my $z = Math::Complex->make( $x, $y);
556 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:
558 my $vector1 = $matrix->col(0)->as_vector;
559 my $vector2 = $matrix->row(1)->as_vector;
566 croak(__PACKAGE__.'::as_vector() - must be a single row or column matrix') unless ($self->cols == 1 or $self->rows == 1);
568 # TODO: there is a faster way to do this
569 return Math::GSL::Vector->new([ $self->as_list ] );
575 Returns a matrix as a flattened Perl list.
577 my @values = $matrix->as_list;
586 for my $row ( 0 .. $self->rows-1) {
587 push @matrix, map { gsl_matrix_get($self->raw, $row, $_) } (0 .. $self->cols-1 );
594 Returns a row matrix of the row you enter.
596 my $matrix = Math::GSL::Matrix->new(3,3);
598 my $matrix_row = $matrix->row(0);
604 my ($self, $row) = @_;
605 croak (__PACKAGE__.'::$matrix->row($row) - invalid $row value')
606 unless (($row < $self->rows) and $row >= 0);
608 my $rowmat = Math::GSL::Matrix->new(1,$self->cols);
610 for my $n (0 .. $self->cols-1) {
611 gsl_matrix_set( $rowmat->raw, 0, $n,
612 gsl_matrix_get($self->raw, $row, $n)
621 Returns a col matrix of the column you enter.
623 my $matrix = Math::GSL::Matrix->new(3,3);
625 my $matrix_col = $matrix->col(0);
631 my ($self, $col) = @_;
632 croak (__PACKAGE__."::\$matrix->col(\$col) - $col not a valid column")
633 unless ($col < $self->cols and $col >= 0);
635 my $colmat = Math::GSL::Matrix->new($self->rows, 1);
637 map { gsl_matrix_set($colmat->raw, $_, 0,
638 gsl_matrix_get($self->raw, $_, $col),
640 } (0 .. $self->rows-1);
647 Sets a the values of a row with the elements of an array.
649 my $matrix = Math::GSL::Matrix->new(3,3);
650 $matrix->set_row(0, [8, 6, 2]);
652 You can also set multiple rows at once with chained calls:
653 my $matrix = Math::GSL::Matrix->new(3,3);
654 $matrix->set_row(0, [8, 6, 2])
655 ->set_row(1, [2, 4, 1]);
661 my ($self, $row, $values) = @_;
662 my $length = $#$values;
663 die __PACKAGE__.'::new($x, $values) - $values must be a nonempty array reference' if $length == -1;
664 die __PACKAGE__.'::set_row($x, $values) - $x must be a valid row number'
665 if ($row < 0 || $row >= $self->rows);
666 die __PACKAGE__.'::set_row($x, $values) - $values must contains the same number of elements as there is columns in the matrix'
667 if($length != $self->cols-1);
668 map { gsl_matrix_set($self->raw, $row, $_, $values->[$_]) } (0..$length);
674 Sets a the values of a column with the elements of an array.
676 my $matrix = Math::GSL::Matrix->new(3,3);
677 $matrix->set_col(0, [8, 6, 2]);
679 You can also set multiple columns at once with chained calls:
680 my $matrix = Math::GSL::Matrix->new(3,3);
681 $matrix->set_col(0, [8, 6, 2])
682 ->set_col(1, [2, 4, 1]);
688 my ($self, $col, $values) = @_;
689 my $length = $#$values;
691 die __PACKAGE__.'::new($x, $values) - $values must be a nonempty array reference' if $length == -1;
692 die __PACKAGE__.'::set_row($x, $values) - $x must be a valid column number'
693 if ($col < 0 || $col >= $self->cols);
694 die __PACKAGE__.'::set_row($x, $values) - $values must contains the same number of elements as there is rowss in the matrix'
695 if($length != $self->rows-1);
697 map { gsl_matrix_set($self->raw, $_, $col, $values->[$_]) } (0..$length);
702 my ($left, $right) = @_;
704 my $lcopy = $left->copy;
706 if ( blessed $right && $right->isa('Math::GSL::Matrix') && blessed $left && $left->isa('Math::GSL::Matrix') ) {
707 if ( $left->cols == $right->cols && $left->rows == $right->rows ) {
708 gsl_matrix_add($lcopy->raw, $right->raw);
710 croak "Math::GSL - addition of matrices must be called with two objects matrices and must have the same number of columns and rows";
713 gsl_matrix_add_constant($lcopy->raw, $right);
719 my ($left, $right) = @_;
721 my $lcopy = $left->copy;
723 if ( blessed $right && $right->isa('Math::GSL::Matrix') && blessed $left && $left->isa('Math::GSL::Matrix') ) {
724 if ( $left->cols == $right->cols && $left->rows == $right->rows ) {
725 gsl_matrix_sub($lcopy->raw, $right->raw);
727 croak "Math::GSL - subtraction of matrices must be called with two objects matrices and must have the same number of columns and rows";
730 gsl_matrix_add_constant($lcopy->raw, $right*-1);
735 sub _multiplication {
736 my ($left,$right) = @_;
737 my $lcopy = $left->copy;
739 if ( blessed $right && $right->isa('Math::GSL::Matrix') ) {
740 if ( $left->cols == $right->cols && $left->rows == $right->rows ) {
741 return _dot_product($left,$right);
742 #gsl_matrix_mul_elements($lcopy->raw, $right->raw);
744 croak "Math::GSL::Matrix - multiplication of elements of matrices must be called with two objects matrices and must have the same number of columns and rows";
747 gsl_matrix_scale($lcopy->raw, $right);
753 my ($left, $right) = @_;
755 return is_similar( [ $left->as_list ], [ $right->as_list ] );
759 my ($left, $right) = @_;
760 return !is_similar( [ $left->as_list ], [ $right->as_list ] );
764 my ($left,$right) = @_;
766 croak "Math::GSL::Matrix - incompatible matrices in multiplication"
767 unless ( blessed $right && $right->isa('Math::GSL::Matrix') and $left->rows == $right->cols );
768 my $C = Math::GSL::Matrix->new($left->rows, $right->cols);
769 gsl_blas_dgemm($CblasNoTrans, $CblasNoTrans, 1, $left->raw, $right->raw, 1, $C->raw);
775 Here is a list of all the functions included in this module :
779 =item C<gsl_matrix_alloc($i, $j)> - Return a gsl_matrix of $i rows and $j columns
781 =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
783 =item C<gsl_matrix_alloc_from_block> -
785 =item C<gsl_matrix_free> -
787 =item C<gsl_matrix_alloc_from_matrix > -
789 =item C<gsl_vector_alloc_row_from_matrix> -
791 =item C<gsl_vector_alloc_col_from_matrix > -
793 =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.
795 =item C<gsl_matrix_row($m , $i)> - Return a vector view of the $i-th row of the matrix $m
797 =item C<gsl_matrix_column($m, $j)> - Return a vector view of the $j-th column of the matrix $m
799 =item C<gsl_matrix_diagonal($m)> - Return a vector view of the diagonal of the vector. The matrix doesn't have to be square.
801 =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.
803 =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.
805 =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.
807 =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.
809 =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.
811 =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.
813 =item C<gsl_matrix_view_vector> -
815 =item C<gsl_matrix_view_vector_with_tda> -
817 =item C<gsl_matrix_const_submatrix> -
819 =item C<gsl_matrix_get($m, $i, $j)> - Return the (i,j)-th element of the matrix $m
821 =item C<gsl_matrix_set($m, $i, $j, $x)> - Set the value of the (i,j)-th element of the matrix $m to $x
823 =item C<gsl_matrix_ptr> -
825 =item C<gsl_matrix_const_ptr> -
827 =item C<gsl_matrix_set_zero($m)> - Set all the elements of the matrix $m to zero
829 =item C<gsl_matrix_set_identity($m)> - Set the elements of the matrix $m to the corresponding elements of the identity matrix
831 =item C<gsl_matrix_set_all($m, $x)> - Set all the elements of the matrix $m to the value $x
833 =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
835 =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
837 =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
839 =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
841 =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.
843 =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.
845 =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 suceeded, 1 otherwise.
847 =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 suceeded, 1 otherwise.
849 =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 suceeded, 1 otherwise.
851 =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.
853 =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.
855 =item C<gsl_matrix_max($m)> - Return the maximum value in the matrix $m
857 =item C<gsl_matrix_min($m)> - Return the minimum value in the matrix $m
859 =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.
861 =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.
863 =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.
865 =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
867 =item C<gsl_matrix_isnull($m)> - Return 1 if all the elements of the matrix $m are zero, 0 otherwise
869 =item C<gsl_matrix_ispos($m)> - Return 1 if all the elements of the matrix $m are strictly positve, 0 otherwise
871 =item C<gsl_matrix_isneg($m)> - Return 1 if all the elements of the matrix $m are strictly negative, 0 otherwise
873 =item C<gsl_matrix_isnonneg($m)> - Return 1 if all the elements of the matrix $m are non-negatuive, 0 otherwise
875 =item C<gsl_matrix_add($a, $b)> - Add the elements of matrix $b to the elements of matrix $a
877 =item C<gsl_matrix_sub($a, $b)> - Subtract the elements of matrix $b from the elements of matrix $a
879 =item C<gsl_matrix_mul_elements($a, $b)> - Multiplie the elements of matrix $a by the elements of matrix $b
881 =item C<gsl_matrix_div_elements($a, $b)> - Divide the elements of matrix $a by the elements of matrix $b
883 =item C<gsl_matrix_scale($a, $x)> - Multiplie the elements of matrix $a by the constant factor $x
885 =item C<gsl_matrix_add_constant($a, $x)> - Add the constant value $x to the elements of the matrix $a
887 =item C<gsl_matrix_add_diagonal($a, $x)> - Add the constant value $x to the elements of the diagonal of the matrix $a
889 =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 lenght of the vector must be of the same as the lenght of the row. The function returns 0 if it succeded, 1 otherwise.
891 =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 lenght of the vector must be of the same as the lenght of the column. The function returns 0 if it succeded, 1 otherwise.
893 =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 lenght of the vector must be of the same as the lenght of the row. The function returns 0 if it succeded, 1 otherwise.
895 =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 lenght of the vector must be of the same as the lenght of the column. The function returns 0 if it succeded, 1 otherwise.
899 These are related to constant views of a matrix.
904 =item C<gsl_matrix_const_row>
906 =item C<gsl_matrix_const_colum>
908 =item C<gsl_matrix_const_diagonal>
910 =item C<gsl_matrix_const_subdiagonal>
912 =item C<gsl_matrix_const_superdiagonal>
914 =item C<gsl_matrix_const_subrow>
916 =item C<gsl_matrix_const_subcolumn>
918 =item C<gsl_matrix_const_view_array>
920 =item C<gsl_matrix_const_view_array_with_tda>
925 The following functions are similar to those above but work with C<char>'s and C<int>'s. We are not quite
926 sure if anyone wants these. Please speak up if you do and/or submit some patches to this documentation, please!
931 =item gsl_matrix_const_view_vector
933 =item gsl_matrix_const_view_vector_with_tda
935 =item gsl_matrix_char_alloc
937 =item gsl_matrix_char_calloc
939 =item gsl_matrix_char_alloc_from_block
941 =item gsl_matrix_char_alloc_from_matrix
943 =item gsl_vector_char_alloc_row_from_matrix
945 =item gsl_vector_char_alloc_col_from_matrix
947 =item gsl_matrix_char_free
949 =item gsl_matrix_char_submatrix
951 =item gsl_matrix_char_row
953 =item gsl_matrix_char_column
955 =item gsl_matrix_char_diagonal
957 =item gsl_matrix_char_subdiagonal
959 =item gsl_matrix_char_superdiagonal
961 =item gsl_matrix_char_subrow
963 =item gsl_matrix_char_subcolumn
965 =item gsl_matrix_char_view_array
967 =item gsl_matrix_char_view_array_with_tda
969 =item gsl_matrix_char_view_vector
971 =item gsl_matrix_char_view_vector_with_tda
973 =item gsl_matrix_char_const_submatrix
975 =item gsl_matrix_char_const_row
977 =item gsl_matrix_char_const_column
979 =item gsl_matrix_char_const_diagonal
981 =item gsl_matrix_char_const_subdiagonal
983 =item gsl_matrix_char_const_superdiagonal
985 =item gsl_matrix_char_const_subrow
987 =item gsl_matrix_char_const_subcolumn
989 =item gsl_matrix_char_const_view_array
991 =item gsl_matrix_char_const_view_array_with_tda
993 =item gsl_matrix_char_const_view_vector
995 =item gsl_matrix_char_const_view_vector_with_tda
997 =item gsl_matrix_char_get
999 =item gsl_matrix_char_set
1001 =item gsl_matrix_char_ptr
1003 =item gsl_matrix_char_const_ptr
1005 =item gsl_matrix_char_set_zero
1007 =item gsl_matrix_char_set_identity
1009 =item gsl_matrix_char_set_all
1011 =item gsl_matrix_char_fread
1013 =item gsl_matrix_char_fwrite
1015 =item gsl_matrix_char_fscanf
1017 =item gsl_matrix_char_fprintf
1019 =item gsl_matrix_char_memcpy
1021 =item gsl_matrix_char_swap
1023 =item gsl_matrix_char_swap_rows
1025 =item gsl_matrix_char_swap_columns
1027 =item gsl_matrix_char_swap_rowcol
1029 =item gsl_matrix_char_transpose
1031 =item gsl_matrix_char_transpose_memcpy
1033 =item gsl_matrix_char_max
1035 =item gsl_matrix_char_min
1037 =item gsl_matrix_char_minmax
1039 =item gsl_matrix_char_max_index
1041 =item gsl_matrix_char_min_index
1043 =item gsl_matrix_char_minmax_index
1045 =item gsl_matrix_char_isnull
1047 =item gsl_matrix_char_ispos
1049 =item gsl_matrix_char_isneg
1051 =item gsl_matrix_char_isnonneg
1053 =item gsl_matrix_char_add
1055 =item gsl_matrix_char_sub
1057 =item gsl_matrix_char_mul_elements
1059 =item gsl_matrix_char_div_elements
1061 =item gsl_matrix_char_scale
1063 =item gsl_matrix_char_add_constant
1065 =item gsl_matrix_char_add_diagonal
1067 =item gsl_matrix_int_alloc
1069 =item gsl_matrix_int_calloc
1071 =item gsl_matrix_int_alloc_from_block
1073 =item gsl_matrix_int_alloc_from_matrix
1075 =item gsl_vector_int_alloc_row_from_matrix
1077 =item gsl_vector_int_alloc_col_from_matrix
1079 =item gsl_matrix_int_free
1081 =item gsl_matrix_int_submatrix
1083 =item gsl_matrix_int_row
1085 =item gsl_matrix_int_column
1087 =item gsl_matrix_int_diagonal
1089 =item gsl_matrix_int_subdiagonal
1091 =item gsl_matrix_int_superdiagonal
1093 =item gsl_matrix_int_subrow
1095 =item gsl_matrix_int_subcolumn
1097 =item gsl_matrix_int_view_array
1099 =item gsl_matrix_int_view_array_with_tda
1101 =item gsl_matrix_int_view_vector
1103 =item gsl_matrix_int_view_vector_with_tda
1105 =item gsl_matrix_int_const_submatrix
1107 =item gsl_matrix_int_const_row
1109 =item gsl_matrix_int_const_column
1111 =item gsl_matrix_int_ptr
1113 =item gsl_matrix_int_const_ptr
1115 =item gsl_matrix_int_set_zero
1117 =item gsl_matrix_int_set_identity
1119 =item gsl_matrix_int_set_all
1121 =item gsl_matrix_int_fread
1123 =item gsl_matrix_int_fwrite
1125 =item gsl_matrix_int_fscanf
1127 =item gsl_matrix_int_fprintf
1129 =item gsl_matrix_int_memcpy
1131 =item gsl_matrix_int_swap
1133 =item gsl_matrix_int_swap_rows
1135 =item gsl_matrix_int_swap_columns
1137 =item gsl_matrix_int_swap_rowcol
1139 =item gsl_matrix_int_transpose
1141 =item gsl_matrix_int_transpose_memcpy
1143 =item gsl_matrix_int_max
1145 =item gsl_matrix_int_min
1147 =item gsl_matrix_int_minmax
1149 =item gsl_matrix_int_max_index
1151 =item gsl_matrix_int_min_index
1153 =item gsl_matrix_int_minmax_index
1155 =item gsl_matrix_int_isnull
1157 =item gsl_matrix_int_ispos
1159 =item gsl_matrix_int_isneg
1161 =item gsl_matrix_int_isnonneg
1163 =item gsl_matrix_int_add
1165 =item gsl_matrix_int_sub
1167 =item gsl_matrix_int_mul_elements
1169 =item gsl_matrix_int_div_elements
1171 =item gsl_matrix_int_scale
1173 =item gsl_matrix_int_add_constant
1175 =item gsl_matrix_int_add_diagonal
1179 You have to add the functions you want to use inside the qw /put_funtion_here /.
1180 You can also write use Math::GSL::Matrix qw/:all/ to use all avaible functions of the module.
1181 Other tags are also avaible, here is a complete list of all tags for this module :
1197 For more informations on the functions, we refer you to the GSL offcial documentation
1198 L<http://www.gnu.org/software/gsl/manual/html_node/>
1200 Tip : search on google: site:http://www.gnu.org/software/gsl/manual/html_node/name_of_the_function_you_want
1204 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>
1206 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.
1208 use Math::GSL::Matrix qw/:all/;
1209 my $m = gsl_matrix_alloc (10,3);
1212 gsl_matrix_set($m, $i, $j, 0.23 + 100*$i + $j);
1216 for my $i (0..99){ # OUT OF RANGE ERROR
1218 print "m($i, $j) = " . gsl_matrix_get ($m, $i, $j) . "\n";
1221 gsl_matrix_free ($m);
1224 use Math::GSL::Matrix qw/:all/;
1226 my $m = gsl_matrix_alloc (100, 100);
1227 my $a = gsl_matrix_alloc (100, 100);
1231 gsl_matrix_set ($m, $i, $j, 0.23 + $i + $j);
1235 The next program shows how to write a matrix to a file.
1237 my $out = gsl_fopen("test.dat", "wb");
1238 gsl_matrix_fwrite ($out, $m);
1241 my $in = gsl_fopen("test.dat", "rb");
1242 gsl_matrix_fread ($in, $a);
1248 $mij = gsl_matrix_get ($m, $i, $j);
1249 $aij = gsl_matrix_get ($a, $i, $j);
1250 $k++ if ($mij != $aij);
1254 gsl_matrix_free($m);
1255 gsl_matrix_free($a);
1257 print "differences = $k (should be zero)\n";
1261 Jonathan Leto <jonathan@leto.net> and Thierry Moisan <thierry.moisan@gmail.com>
1263 =head1 COPYRIGHT AND LICENSE
1265 Copyright (C) 2008-2009 Jonathan Leto and Thierry Moisan
1267 This program is free software; you can redistribute it and/or modify it
1268 under the same terms as Perl itself.