pod/Vector.pod

   1 %perlcode %{
   2 use Scalar::Util 'blessed';
   3 use Data::Dumper;
   4 use Carp qw/croak/;
   5 use Math::GSL::Errno qw/:all/;
   6 use Math::GSL::BLAS qw/gsl_blas_ddot/;
   7 use overload
   8     '*'      => \&_multiplication,
   9     '+'      => \&_addition,
  10     '-'      => \&_subtract,
  11     fallback => 1,
  12 ;
  13
  14 @EXPORT_all  = qw/fopen fclose
  15                  gsl_vector_alloc gsl_vector_calloc gsl_vector_alloc_from_block gsl_vector_alloc_from_vector
  16                  gsl_vector_free gsl_vector_view_array gsl_vector_const_view_array gsl_vector_view_array_with_stride
  17                  gsl_vector_const_view_array_with_stride gsl_vector_subvector gsl_vector_subvector_wi gsl_vector_subvector_with_stride
  18                  gsl_vector_const_subvec gsl_vector_const_subvec gsl_vector_get gsl_vector_set
  19                  gsl_vector_ptr gsl_vector_const_ptr gsl_vector_set_zero gsl_vector_set_all
  20                  gsl_vector_set_basis gsl_vector_fread gsl_vector_fwrite gsl_vector_fscanf
  21                  gsl_vector_fprintf gsl_vector_memcpy gsl_vector_reverse gsl_vector_swap
  22                  gsl_vector_swap_elements gsl_vector_max gsl_vector_min gsl_vector_minmax
  23                  gsl_vector_max_index gsl_vector_min_index gsl_vector_minmax_index
  24                  gsl_vector_add gsl_vector_sub gsl_vector_mul gsl_vector_div
  25                  gsl_vector_scale gsl_vector_add_constant gsl_vector_isnull
  26                  gsl_vector_ispos gsl_vector_isneg gsl_vector_isnonneg
  27                  gsl_vector_float_alloc gsl_vector_float_calloc gsl_vector_float_alloc_from_block
  28                  gsl_vector_float_alloc_from_vector gsl_vector_float_free gsl_vector_float_view_array
  29                  gsl_vector_float_view_array_with_stride gsl_vector_float_const_view_array gsl_vector_float_const_view_array_with_stride
  30                  gsl_vector_float_subvector gsl_vector_float_subvector_with_stride gsl_vector_float_const_subvector
  31                  gsl_vector_float_const_subvector_with_stride gsl_vector_float_get gsl_vector_float_set gsl_vector_float_ptr
  32                  gsl_vector_float_const_ptr gsl_vector_float_set_zero gsl_vector_float_set_all gsl_vector_float_set_basis
  33                  gsl_vector_float_fread gsl_vector_float_fwrite gsl_vector_float_fscanf gsl_vector_float_fprintf
  34                  gsl_vector_float_memcpy gsl_vector_float_reverse gsl_vector_float_swap gsl_vector_float_swap_elements
  35                  gsl_vector_float_max gsl_vector_float_min gsl_vector_float_minmax gsl_vector_float_max_index gsl_vector_float_min_index
  36                  gsl_vector_float_minmax_index gsl_vector_float_add gsl_vector_float_sub gsl_vector_float_mul gsl_vector_float_div gsl_vector_float_scale
  37                  gsl_vector_float_add_constant gsl_vector_float_isnull gsl_vector_float_ispos gsl_vector_float_isneg gsl_vector_float_isnonneg
  38                  gsl_vector_complex_alloc gsl_vector_complex_calloc gsl_vector_complex_alloc_from_block gsl_vector_complex_alloc_from_vector
  39                  gsl_vector_complex_free gsl_vector_complex_view_array gsl_vector_complex_view_array_with_stride gsl_vector_complex_const_view_array
  40                  gsl_vector_complex_const_view_array_with_stride gsl_vector_complex_subvector gsl_vector_complex_subvector_with_stride
  41                  gsl_vector_complex_const_subvector gsl_vector_complex_const_subvector_with_stride gsl_vector_complex_real gsl_vector_complex_imag
  42                  gsl_vector_complex_const_real gsl_vector_complex_const_imag gsl_vector_complex_get gsl_vector_complex_set
  43                  gsl_vector_complex_ptr gsl_vector_complex_const_ptr gsl_vector_complex_set_zero gsl_vector_complex_set_all
  44                  gsl_vector_complex_set_basis gsl_vector_complex_fread gsl_vector_complex_fwrite gsl_vector_complex_fscanf
  45                  gsl_vector_complex_fprintf gsl_vector_complex_memcpy gsl_vector_complex_reverse gsl_vector_complex_swap
  46                  gsl_vector_complex_swap_elements gsl_vector_complex_isnull gsl_vector_complex_ispos gsl_vector_complex_isneg
  47 /;
  48 @EXPORT_file =qw/ fopen fclose/;
  49 @EXPORT_OK = (@EXPORT_all, @EXPORT_file);
  50 %EXPORT_TAGS = ( file => \@EXPORT_file, all => \@EXPORT_all );
  51
  52 =head1 NAME
  53
  54 Math::GSL::Vector - Functions concerning vectors
  55
  56 =head1 SYNOPSIS
  57
  58     use Math::GSL::Vector qw/:all/;
  59     my $vec1 = Math::GSL::Vector->new([1, 7, 94, 15 ]);
  60     my $vec2 = $vec1 * 5;
  61     my $vec3 = Math::GSL::Vector>new(10);   # 10 element zero vector
  62     my $vec4 = $vec1 + $vec2;
  63
  64     # set the element at index 1 to 9
  65     # and the element at index 3 to 8
  66     $vec3->set([ 1, 3 ], [ 9, 8 ]);
  67
  68     my @vec = $vec2->as_list;               # return elements as Perl list
  69
  70     my $dot_product = $vec1 * $vec2;
  71     my $length      = $vec2->length;
  72     my $first       = $vec1->get(0);
  73
  74
  75 =cut
  76
  77 =head1 Objected Oriented Interface to GSL Math::GSL::Vector
  78
  79 =head2 Math::GSL::Vector->new()
  80
  81 Creates a new Vector of the given size.
  82
  83     my $vector = Math::GSL::Matrix->new(3);
  84
  85 You can also create and set directly the values of the vector like this :
  86
  87    my $vector = Math::GSL::Vector->new([2,4,1]);
  88
  89 =cut
  90 sub new {
  91     my ($class, $values) = @_;
  92     my $length  = $#$values;
  93     my $this = {};
  94     my $vector;
  95     if ( ref $values eq 'ARRAY' ){
  96         die __PACKAGE__.'::new($x) - $x must be a nonempty array reference' if $length == -1;
  97         $vector  = gsl_vector_alloc($length+1);
  98         map { gsl_vector_set($vector, $_, $values->[$_] ) }  (0 .. $length);
  99         $this->{_length} = $length+1;
 100     } elsif ( (int($values) == $values) && ($values > 0)) {
 101         $vector  = gsl_vector_alloc($values);
 102         gsl_vector_set_zero($vector);
 103         $this->{_length} = $values;
 104     } else {
 105         die __PACKAGE__.'::new($x) - $x must be an int or array reference';
 106     }
 107     $this->{_vector} = $vector;
 108     bless $this, $class;
 109 }
 110 =head2 raw()
 111
 112 Get the underlying GSL vector object created by SWIG, useful for using gsl_vector_* functions which do not have an OO counterpart.
 113
 114     my $vector    = Math::GSL::Vector->new(3);
 115     my $gsl_vector = $vector->raw;
 116     my $stuff      = gsl_vector_get($gsl_vector, 1);
 117
 118 =cut
 119
 120 sub raw {
 121     my $self = shift;
 122     return $self->{_vector};
 123 }
 124
 125 =head2 min()
 126
 127 Returns the minimum value contained in the vector.
 128
 129    my $vector = Math::GSL::Vector->new([2,4,1]);
 130    my $minimum = $vector->min;
 131
 132 =cut
 133
 134 sub min {
 135     my $self=shift;
 136     return gsl_vector_min($self->raw);
 137 }
 138
 139 =head2 max()
 140
 141 Returns the minimum value contained in the vector.
 142
 143    my $vector = Math::GSL::Vector->new([2,4,1]);
 144    my $maximum = $vector->max;
 145
 146 =cut
 147
 148 sub max {
 149     my $self=shift;
 150     return gsl_vector_max($self->raw);
 151 }
 152
 153 =head2 length()
 154
 155 Returns the number of elements contained in the vector.
 156
 157    my $vector = Math::GSL::Vector->new([2,4,1]);
 158    my $length = $vector->length;
 159
 160 =cut
 161
 162 sub length { my $self=shift; $self->{_length} }
 163
 164 =head2  as_list()
 165
 166 Gets the content of a Math::GSL::Vector object as a Perl list.
 167
 168     my $vector = Math::GSL::Vector->new(3);
 169     ...
 170     my @values = $vector->as_list;
 171 =cut
 172
 173 sub as_list {
 174     my $self=shift;
 175     $self->get( [ 0 .. $self->length - 1  ] );
 176 }
 177
 178 =head2  get()
 179
 180 Gets the value of an of a Math::GSL::Vector object.
 181
 182     my $vector = Math::GSL::Vector->new(3);
 183     ...
 184     my @values = $vector->get(2);
 185
 186 You can also enter an array of indices to receive their corresponding values:
 187
 188     my $vector = Math::GSL::Vector->new(3);
 189     ...
 190     my @values = $vector->get([0,2]);
 191
 192 =cut
 193
 194 sub get {
 195     my ($self, $indices) = @_;
 196     return  map {  gsl_vector_get($self->raw, $_ ) } @$indices ;
 197 }
 198
 199 =head2  set()
 200
 201 Sets values of an of a Math::GSL::Vector object.
 202
 203     my $vector = Math::GSL::Vector->new(3);
 204     $vector->set([1,2], [8,23]);
 205
 206 This sets the second and third value to 8 and 23.
 207
 208 =cut
 209
 210 sub set {
 211     my ($self, $indices, $values) = @_;
 212     die (__PACKAGE__.'::set($indices, $values) - $indices and $values must be array references of the same length')
 213         unless ( ref $indices eq 'ARRAY' && ref $values eq 'ARRAY' &&  $#$indices == $#$values );
 214     eval {
 215         map {  gsl_vector_set($self->{_vector}, $indices->[$_], $values->[$_] ) } (0..$#$indices);
 216     };
 217     return;
 218 }
 219
 220 =head2 copy()
 221
 222 Returns a copy of the vector, which has the same length and values but resides at a different location in memory.
 223
 224     my $vector = Math::GSL::Vector->new([10 .. 20]);
 225     my $copy   = $vector->copy;
 226
 227 =cut
 228
 229
 230 sub copy {
 231     my $self = shift;
 232     my $copy = Math::GSL::Vector->new( $self->length );
 233     if ( gsl_vector_memcpy($copy->raw, $self->raw) != $GSL_SUCCESS ) {
 234         croak "Math::GSL - error copying memory, aborting";
 235     }
 236     return $copy;
 237 }
 238
 239 sub _multiplication {
 240     my ($left,$right) = @_;
 241     my $lcopy = $left->copy;
 242
 243     if ( blessed $right && $right->isa('Math::GSL::Vector') ) {
 244         return $lcopy->dot_product($right);
 245     } else {
 246         gsl_vector_scale($lcopy->raw, $right);
 247     }
 248     return $lcopy;
 249 }
 250
 251 sub _subtract {
 252     my ($left, $right, $flip) = @_;
 253
 254     if ($flip) {
 255         my $lcopy = $left->copy;
 256         gsl_vector_scale($lcopy->raw, -1 );
 257         gsl_vector_add_constant($lcopy->raw, $right);
 258         return $lcopy;
 259     } else {
 260         return _addition($left, -1.0*$right);
 261     }
 262 }
 263
 264 sub _addition {
 265     my ($left, $right, $flip) = @_;
 266
 267     my $lcopy = $left->copy;
 268
 269     if ( blessed $right && $right->isa('Math::GSL::Vector') && blessed $left && $left->isa('Math::GSL::Vector') ) {
 270         if ( $left->length == $right->length ) {
 271             gsl_vector_add($lcopy->raw, $right->raw);
 272         } else {
 273             croak "Math::GSL - addition of vectors must be called with two objects vectors and must have the same length";
 274         }
 275     } else {
 276         gsl_vector_add_constant($lcopy->raw, $right);
 277     }
 278     return $lcopy;
 279 }
 280
 281 sub dot_product_pp {
 282     my ($left,$right) = @_;
 283     my $sum=0;
 284     if ( blessed $right && $right->isa('Math::GSL::Vector') &&
 285          $left->length == $right->length ) {
 286          my @l = $left->as_list;
 287          my @r = $right->as_list;
 288          map { $sum += $l[$_] * $r[$_] } (0..$#l);
 289         return $sum;
 290     } else {
 291         croak "dot_product() must be called with two vectors";
 292     }
 293 }
 294
 295 sub dot_product {
 296     my ($left,$right) = @_;
 297
 298     my ($status, $product) = gsl_blas_ddot($left->raw,$right->raw);
 299     croak sprintf "Math::GSL::dot_product - %s", gsl_strerror($status) if ($status != $GSL_SUCCESS);
 300     return $product;
 301 }
 302
 303 =head1 DESCRIPTION
 304
 305 Here is a list of all the functions included in this module :
 306
 307 =over 1
 308
 309 =item C<gsl_vector_alloc($x)> - create a vector of size $x
 310
 311 =item C<gsl_vector_calloc($x)> - create a vector of size $x and initializes all the elements of the vector to zero
 312
 313 =item C<gsl_vector_alloc_from_block>
 314
 315 =item C<gsl_vector_alloc_from_vector>
 316
 317 =item C<gsl_vector_free($v)> - free a previously allocated vector $v
 318
 319 =item C<gsl_vector_view_array($base, $n)> - This function returns a vector view of an array reference $base. The start of the new vector is given by $base and has $n elements. Mathematically, the i-th element of the new vector v' is given by, v'(i) = $base->[i] where the index i runs from 0 to $n-1. The array containing the elements of v is not owned by the new vector view. When the view goes out of scope the original array will continue to exist. The original memory can only be deallocated by freeing the original pointer base. Of course, the original array should not be deallocated while the view is still in use.
 320
 321 =item C<gsl_vector_const_view_array($base, $n)> - This function is equivalent to gsl_vector_view_array but can be used for arrays which are declared const.
 322
 323 =item C<gsl_vector_view_array_with_stride($base, $stride, $n)> - This function returns a vector view of an array reference $base with an additional $stride argument. The subvector is formed in the same way as for gsl_vector_view_array but the new vector has $n elements with a step-size of $stride from one element to the next in the original array. Mathematically, the i-th element of the new vector v' is given by, v'(i) = $base->[i*$stride] where the index i runs from 0 to $n-1. Note that the view gives direct access to the underlying elements of the original array. A vector view $view can be passed to any subroutine which takes a vector argument just as a directly allocated vector would be, using $view->{vector}.
 324
 325 =item C<gsl_vector_const_view_array_with_stride($base, $stride, $n)> - This function is equivalent to gsl_vector_view_array_with_stride but can be used for arrays which are declared const.
 326
 327 =item C<gsl_vector_subvector($v, $offset, $n)> - return a vector_view type which contains a subvector of $v, with a size of $size, starting from the $offset position
 328
 329 =item C<gsl_vector_subvector_with_stride($v, $offset, $stride, $size)> - return a vector_view type which contains a subvector of $v, with a size of $size, starting from the $offset position and with a $stride step between each element of $v
 330
 331 =item C<gsl_vector_const_subvector>
 332
 333 =item C<gsl_vector_get($v, $i)> - return the $i-th element of a vector $v
 334
 335 =item C<gsl_vector_set($v, $i, $x)> - return the vector $v with his $i-th element set to $x
 336
 337 =item C<gsl_vector_ptr>
 338
 339 =item C<gsl_vector_const_ptr>
 340
 341 =item C<gsl_vector_set_zero($v)> - set all the elements of $v to 0
 342
 343 =item C<gsl_vector_set_all($v, $x)> - set all the elements of $v to $x
 344
 345 =item C<gsl_vector_set_basis($v, $i)> - set all the elements of $v to 0 except for the $i-th element which is set to 1 and return 0 if the operation succeded, 1 otherwise.
 346
 347 =item C<gsl_vector_fread($file, $v)> - This function reads into the vector $v from the open stream $file opened with gsl_fopen function from the Math::GSL module in binary format. The vector $v must be preallocated with the correct length since the function uses the size of $v to determine how many bytes to read. The return value is 0 for success and 1 if there was a problem reading from the file.
 348
 349 =item C<gsl_vector_fwrite($file, $v)> - This function writes the elements of the vector $v to the stream $file opened with gsl_fopen function from the Math::GSL module in binary format. The return value is 0 for success and 1 if there was a problem writing to the file. Since the data is written in the native binary format it may not be portable between different architectures.
 350
 351 =item C<gsl_vector_fscanf($file, $v)> This function reads formatted data from the stream $file opened with gsl_fopen function from the Math::GSL module into the vector $v. The vector $v must be preallocated with the correct length since the function uses the size of $v to determine how many numbers to read. The function returns 0 for success and 1 if there was a problem reading from the file.
 352
 353 =item C<gsl_vector_fprintf($file, $v, $format)> -This function writes the elements of the vector $v line-by-line to the stream $file opened with gsl_fopen function from the Math::GSL module using the format specifier $format, which should be one of the "%g", "%e" or "%f" formats for floating point numbers and "%d" for integers. The function returns 0 for success and 1 if there was a problem writing to the file.
 354
 355 =item C<gsl_vector_memcpy($dest, $src)> - This function copies the elements of the vector $src into the vector $dest and return 0 if the opertaion succeded, 1 otherwise. The two vectors must have the same length.
 356
 357 =item C<gsl_vector_reverse($v)> - reverse the order of the elements of the vector $v and return 0 if the opertaion succeded, 1 otherwise
 358
 359 =item C<gsl_vector_swap($v, $v2)> - swap the values of the vectors $v and $v2 and return 0 if the opertaion succeded, 1 otherwise
 360
 361 =item C<gsl_vector_swap_elements($v, $i, $j)> - permute the elements at position $i and $j in the vector $v and return 0 if the operation succeded, 1 otherwise.
 362
 363 =item C<gsl_vector_max($v)> - return the maximum value in the vector $v
 364
 365 =item C<gsl_vector_min($v)> - return the minimum value in the vector $v
 366
 367 =item C<gsl_vector_minmax($v)> - return two values, the first is the minimum value in the vector $v and the second is the maximum value.
 368
 369 =item C<gsl_vector_max_index($v)> - return the position of the maximum value in the vector $v
 370
 371 =item C<gsl_vector_min_index($v)> - return the position of the minimum value in the vector $v
 372
 373 =item C<gsl_vector_minmax_index> -  return two values, the first is the position of the minimum value in the vector $v and the second is the position of the maximum value.
 374
 375 =item C<gsl_vector_add($v, $v2)> - add the elements of $v2 to the elements of $v, the two vectors must have the same lenght and return 0 if the operation succeded, 1 otherwise.
 376
 377 =item C<gsl_vector_sub($v, $v2)> - substract the elements of $v2 from the elements of $v, the two vectors must have the same lenght and return 0 if the operation succeded, 1 otherwise.
 378
 379 =item C<gsl_vector_mul($v, $v2)> - multiply the elements of $v by the elements of $v2, the two vectors must have the same lenght and return 0 if the operation succeded, 1 otherwise.
 380
 381 =item C<gsl_vector_div($v, $v2)> - divides the elements of $v by the elements of $v2, the two vectors must have the same lenght and return 0 if the operation succeded, 1 otherwise.
 382
 383 =item C<gsl_vector_scale($v, $x)> - multiplty the elements of the vector $v by a constant $x and return 0 if the operation succeded, 1 otherwise.
 384
 385 =item C<gsl_vector_add_constant($v, $x)> - add a constant $x to the elements of the vector $v and return 0 if the operation succeded, 1 otherwise.
 386
 387 =item C<gsl_vector_isnull($v)> - verify if all the elements of the vector $v are null, return 0 if it's the case, 1 otherwise.
 388
 389 =item C<gsl_vector_ispos($v)> - verify if all the elements of the vector $v are positive, return 0 if it's the case, 1 otherwise.
 390
 391 =item C<gsl_vector_isneg($v)> - verify if all the elements of the vector $v are negative, return 0 if it's the case, 1 otherwise.
 392
 393 =item C<gsl_vector_isnonneg($v)> - verify if all the elements the vector $v are not negative, return 0 if it's the case, 1 otherwise.
 394
 395 =back
 396
 397 You have to add the functions you want to use inside the qw /put_funtion_here / with spaces between each function. You can also write use Math::GSL::Complex qw/:all/ to use all avaible functions of the module.
 398
 399
 400 Precision on the vector_view type : every modification you'll make on a vector_view will also modify the original vector.
 401 For example, the following code will zero the even elements of the vector $v of length $size, while leaving the odd elements untouched :
 402
 403 =over 1
 404
 405 =item C<$v_even= gsl_vector_subvector_with_stride ($v, 0, 2, $size/2);>
 406
 407 =item C<gsl_vector_set_zero ($v_even-E<gt>{vector});>
 408
 409 =back
 410
 411 For more informations on the functions, we refer you to the GSL offcial documentation:
 412 L<http://www.gnu.org/software/gsl/manual/html_node/>
 413
 414 Tip : search on google: site:http://www.gnu.org/software/gsl/manual/html_node/ name_of_the_function_you_want
 415
 416 =head1 EXAMPLES
 417
 418 Here is an example using both interfaces.
 419
 420  use Math::GSL::Vector qw/:all/;
 421
 422  print "We'll create this vector : [0,1,4,9,16] \n";
 423  my $vector = Math::GSL::Vector->new([0,1,4,9,16]);
 424  my ($min, $max) = gsl_vector_minmax_index($vector->raw);
 425
 426  print "We then check the index value of the maximum and minimum values of the vector. \n";
 427  print "The index of the maximum should be 4 and we received $max \n";
 428  print "The index of the minimum should be 0 and we received $min \n";
 429  print "We'll then swap the first and the third elements of the vector \n";
 430
 431  gsl_vector_swap_elements($vector->raw, 0, 3);
 432  my @got = $vector->as_list;
 433  print "The vector should now be like this : [9,1,4,0,16] \n";
 434  print "and we received : [ @got ]\n";
 435
 436 =head1 AUTHORS
 437
 438 Jonathan Leto <jonathan@leto.net> and Thierry Moisan <thierry.moisan@gmail.com>
 439
 440 =head1 COPYRIGHT AND LICENSE
 441
 442 Copyright (C) 2008 Jonathan Leto and Thierry Moisan
 443
 444 This program is free software; you can redistribute it and/or modify it
 445 under the same terms as Perl itself.
 446
 447 =cut
 448 %}
 449