Skip site navigation (1)Skip section navigation (2)

FreeBSD Manual Pages

  
 
  

home | help
Math::GSL::Matrix(3)  User Contributed Perl Documentation Math::GSL::Matrix(3)

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

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	substract 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

DESCRIPTION
       This module is part of the 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.

OBJECT ORIENTED	API
   Constructor
       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 $m is NOT copied. The new object will just refer to the old
       gsl_matrix.

   Getters
       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);

       dim()

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

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

       Basically a shortcut to "rows" and "cols" methods.

       rows()

       Returns the number of rows in the matrix.

	   my $rows = $matrix->rows;

       cols()

       Returns the number of columns in	the matrix.

	   my $cols = $matrix->cols;

       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;

       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 "col()" or	"row()"	into a vector, like so:

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

       get_elem()

       Returns an element of a matrix.

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

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

       row()

       Returns a row matrix of the row you enter.

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

       col()

       Returns a col matrix of the column you enter.

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

       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();

       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();

   Setters
       identity()

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

	   my $I = $matrix->identity;

       zero()

       Set a matrix to the zero	matrix.

	   $matrix->zero;

       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);

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

       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]);
	   ...

       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]);
	   ...

   Utility Functions
       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;

       is_square()

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

       det()

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

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

       write

       Saves 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
       "gsl_fopen", "gsl_matrix_fwrite"	and "gsl_fclose" manually. Future
       versions	might include helper for this task.

       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.

       read

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

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

       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();

       inverse()

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

	   my $inverse = $matrix->inverse;

       transpose()

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

	   my $transposed = $matrix->transpose;

       eigenvalues()

       Computes	the a matrix eigen values.

       eigenpair()

       ieach()

       Applied a function to each element of a matrix.

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

       This method changes the original	matrix.	See "each" for a non
       destructive method.

       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 new matrix.

       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 new
       matrix.

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

       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
       new matrix.

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

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

       "gsl_matrix_alloc($i, $j)" - Return a gsl_matrix	of $i rows and $j
       columns
       "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
       "gsl_matrix_alloc_from_block" -
       "gsl_matrix_free" -
       "gsl_matrix_alloc_from_matrix " -
       "gsl_vector_alloc_row_from_matrix" -
       "gsl_vector_alloc_col_from_matrix " -
       "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.
       "gsl_matrix_row($m , $i)" - Return a vector view	of the $i-th row of
       the matrix $m
       "gsl_matrix_column($m, $j)" - Return a vector view of the $j-th column
       of the matrix $m
       "gsl_matrix_diagonal($m)" - Return a vector view	of the diagonal	of the
       vector. The matrix doesn't have to be square.
       "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.
       "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.
       "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.
       "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.
       "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.
       "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.
       "gsl_matrix_view_vector"	-
       "gsl_matrix_view_vector_with_tda" -
       "gsl_matrix_const_submatrix" -
       "gsl_matrix_get($m, $i, $j)" - Return the (i,j)-th element of the
       matrix $m
       "gsl_matrix_set($m, $i, $j, $x)"	- Set the value	of the (i,j)-th
       element of the matrix $m	to $x
       "gsl_matrix_ptr"	-
       "gsl_matrix_const_ptr" -
       "gsl_matrix_set_zero($m)" - Set all the elements	of the matrix $m to
       zero
       "gsl_matrix_set_identity($m)" - Set the elements	of the matrix $m to
       the corresponding elements of the identity matrix
       "gsl_matrix_set_all($m, $x)" - Set all the elements of the matrix $m to
       the value $x
       "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
       "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
       "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
       "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
       "gsl_matrix_memcpy($dest, $src)"	- Copy the elements of the matrix $src
       to the matrix $dest. The	two matrices must have the same	size.
       "gsl_matrix_swap($m1, $m2)" - Exchange the elements of the matrices $m1
       and $m2 by copying. The two matrices must have the same size.
       "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.
       "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.
       "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.
       "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.
       "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.
       "gsl_matrix_max($m)" - Return the maximum value in the matrix $m
       "gsl_matrix_min($m)" - Return the minimum value in the matrix $m
       "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.
       "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.
       "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.
       "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
       "gsl_matrix_isnull($m)" - Return	1 if all the elements of the matrix $m
       are zero, 0 otherwise
       "gsl_matrix_ispos($m)" -	Return 1 if all	the elements of	the matrix $m
       are strictly positve, 0 otherwise
       "gsl_matrix_isneg($m)" -	Return 1 if all	the elements of	the matrix $m
       are strictly negative, 0	otherwise
       "gsl_matrix_isnonneg($m)" - Return 1 if all the elements	of the matrix
       $m are non-negatuive, 0 otherwise
       "gsl_matrix_add($a, $b)"	- Add the elements of matrix $b	to the
       elements	of matrix $a
       "gsl_matrix_sub($a, $b)"	- Subtract the elements	of matrix $b from the
       elements	of matrix $a
       "gsl_matrix_mul_elements($a, $b)" - Multiplie the elements of matrix $a
       by the elements of matrix $b
       "gsl_matrix_div_elements($a, $b)" - Divide the elements of matrix $a by
       the elements of matrix $b
       "gsl_matrix_scale($a, $x)" - Multiplie the elements of matrix $a	by the
       constant	factor $x
       "gsl_matrix_add_constant($a, $x)" - Add the constant value $x to	the
       elements	of the matrix $a
       "gsl_matrix_add_diagonal($a, $x)" - Add the constant value $x to	the
       elements	of the diagonal	of the matrix $a
       "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.
       "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.
       "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.
       "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.

       These are related to constant views of a	matrix.

       "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"

       The following functions are similar to those above but work with
       "char"'s	and "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!

       gsl_matrix_const_view_vector
       gsl_matrix_const_view_vector_with_tda
       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_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

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

       "all"
       "int"
       "double"
       "char"
       "complex"

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

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";

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

COPYRIGHT AND LICENSE
       Copyright (C) 2008-2014 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.

perl v5.24.1			  2017-07-02		  Math::GSL::Matrix(3)

NAME | SYNOPSIS | DESCRIPTION | OBJECT ORIENTED API | GSL FUNCTION INTERFACE | EXAMPLES | AUTHORS | COPYRIGHT AND LICENSE

Want to link to this manual page? Use this URL:
<https://www.freebsd.org/cgi/man.cgi?query=Math::GSL::Matrix&sektion=3&manpath=FreeBSD+12.1-RELEASE+and+Ports>

home | help