lib/Math/GSL.pm

   1 package Math::GSL;
   2 use strict;
   3 use warnings;
   4 use Math::GSL::Machine qw/:all/;
   5 use Math::GSL::Const qw/:all/;
   6 use Config;
   7 use Data::Dumper;
   8 use Test::More;
   9 use Scalar::Util qw/looks_like_number/;
  10 require DynaLoader;
  11 require Exporter;
  12 our @ISA = qw(Exporter DynaLoader);
  13 our @EXPORT = qw();
  14 our @EXPORT_OK = qw( is_similar $GSL_MODE_DEFAULT $GSL_PREC_DOUBLE $GSL_PREC_SINGLE $GSL_PREC_APPROX);
  15 our %EXPORT_TAGS = ( all => [ @EXPORT_OK ] );
  16
  17 our ($GSL_PREC_DOUBLE, $GSL_PREC_SINGLE, $GSL_PREC_APPROX ) = 0..2;
  18 our $GSL_MODE_DEFAULT = $GSL_PREC_DOUBLE;
  19
  20 use constant MAX_DOUBLE => 1.7976931348623157e+308;
  21 use constant MIN_DOUBLE => 2.2250738585072014e-308;
  22 use constant MAX_FLOAT  => 3.40282347e+38;
  23 use constant MIN_FLOAT  => 1.175494351e-38;
  24 our $VERSION = 0.043;
  25
  26 =head1 NAME
  27
  28 Math::GSL - Perl interface to the  GNU Scientific Library (GSL) using SWIG
  29
  30 =head1 VERSION
  31
  32 Version 0.43
  33
  34 =cut
  35
  36 =head1 SYNOPSIS
  37
  38     use Math::GSL::RNG qw/:all/;
  39
  40
  41 =head1 AUTHOR
  42
  43 Jonathan Leto, C<< <jonathan@leto.net> >> and Thierry Moisan C<< <thierry.moisan@gmail.com> >>
  44
  45 =head1 BUGS
  46
  47 Please report any bugs or feature requests to the authors directly.
  48
  49
  50 =head1 SUPPORT
  51
  52 You can find documentation for this module with the perldoc command.
  53
  54     perldoc Math::GSL
  55
  56 or online at L<http://leto.net/code/Math-GSL/>
  57
  58 =over 4
  59
  60 =item * AnnoCPAN: Annotated CPAN documentation
  61
  62 L<http://annocpan.org/dist/Math::GSL>
  63
  64 =item * CPAN Ratings
  65
  66 L<http://cpanratings.perl.org/d/Math::GSL>
  67
  68 =item * Search CPAN
  69
  70 L<http://search.cpan.org/dist/Math::GSL>
  71
  72 =back
  73
  74
  75 =head1 ACKNOWLEDGEMENTS
  76
  77
  78 =head1 COPYRIGHT & LICENSE
  79
  80 Copyright 2008 Jonathan Leto, Thierry Moisan all rights reserved.
  81
  82 This program is free software; you can redistribute it and/or modify it
  83 under the same terms as Perl itself.
  84
  85 =cut
  86
  87 sub new
  88 {
  89     my ($self,$args) = @_;
  90     my $class = ref $self || $self || 'Math::GSL';
  91     my $this = { };
  92     bless $this, $class;
  93 }
  94
  95 sub subsystems
  96 {
  97     return qw/
  98         BLAS         Diff         Machine       Permute   Statistics
  99         Block        Eigen        Matrix        Poly      Sum
 100         BSpline      Errno        PowInt        Sys
 101         CBLAS        FFT          Min           IEEEUtils
 102         CDF          Fit          Mode          QRNG      Types
 103         Chebyshev    Function     Monte         RNG       Vector
 104         Heapsort     Multifit     Randist
 105         Combination  Histogram    Multimin      Roots     Wavelet
 106         Complex      Histogram2d  Multiroots    SF        Wavelet2D
 107         Const        Siman
 108         DFT          Integration  NTuple        Sort
 109         DHT          Interp       ODEIV
 110         Deriv        Linalg       Permutation   Spline
 111     /;
 112 }
 113
 114 sub verify_results
 115 {
 116     my ($self,$results,$class) = @_;
 117     my ($x,$val);
 118
 119     while (my($k,$v)=each %$results){
 120         my $eps = 2048*$Math::GSL::Machine::GSL_DBL_EPSILON; # TOL3
 121
 122         defined $class ? ( $x = eval qq{${class}::$k} )
 123                        : ( $x = eval $k);
 124
 125         print $@ if $@;
 126         print "got $x for $k\n" if defined $ENV{DEBUG};
 127
 128         if (ref $v eq 'ARRAY'){
 129             ($val, $eps)   = @$v;
 130         } else {
 131              $val = $v;
 132         }
 133         if (!defined $x ){
 134             ok(0, qq{'$k' died} );
 135         } elsif ($x =~ /nan|inf/i){
 136                 ok( $val eq $x, "'$val'?='$x'" );
 137         } else {
 138             my $res = abs($x-$val);
 139             $@ ? ok(0)
 140             : ok( $res <= $eps, "$k ?= $x,\n+- $res, tol=$eps" );
 141         }
 142     }
 143 }
 144
 145 sub is_similar {
 146     my ($x,$y, $eps) = @_;
 147     $eps ||= 1e-8;
 148     if (ref $x eq 'ARRAY' && ref $y eq 'ARRAY') {
 149         if ( $#$x != $#$y ){
 150             warn "is_similar(): argument of different length!";
 151             return 0;
 152         } else {
 153             map {
 154                     my $delta = abs($x->[$_] - $y->[$_]);
 155                     if($delta > $eps){
 156                         warn "\n\tElements start differing at index $_, delta = $delta\n";
 157                         warn qq{\t\t\$x->[$_] = } . $x->[$_] . "\n";
 158                         warn qq{\t\t\$y->[$_] = } . $y->[$_] . "\n";
 159                         return 0;
 160                     }
 161                 } (0..$#$x);
 162             return 1;
 163         }
 164     } else {
 165         abs($x-$y) <= $eps ? return 1 : return 0;
 166     }
 167 }
 168
 169 sub is_valid_double
 170 {
 171     my $x=shift;
 172     return 0 unless ( defined $x && looks_like_number($x) );
 173
 174     return 1 if ($x == 0);
 175
 176     $x = abs $x;
 177     (
 178       $x > MIN_DOUBLE &&
 179       $x < MAX_DOUBLE
 180     ) ? 1 : 0;
 181 }
 182 sub is_valid_float
 183 {
 184     my $x=shift;
 185     return 0 unless ( defined $x && looks_like_number($x) );
 186
 187     return 1 if ($x == 0);
 188
 189     $x = abs $x ;
 190     (
 191       $x > MIN_FLOAT &&
 192       $x < MAX_FLOAT
 193     ) ? 1 : 0;
 194 }
 195
 196 sub _has_quads          { $Config{use64bitint} eq 'define' || ($Config{longsize} >= 8) }
 197 sub _has_long_doubles                 { $Config{d_longdbl}     eq 'define'             }
 198 sub _has_long_doubles_as_default      { $Config{uselongdouble} eq 'define'             }
 199 sub _has_long_doubles_same_as_doubles { $Config{doublesize}    == $Config{longdblsize} }
 200
 201 sub _assert_dies($;$)
 202 {
 203     my ($code,$msg) = @_;
 204     my $status = eval { &$code };
 205     print "status=||$status||\n\$\?=$?\n\$\!=$!\n" if 0;
 206     $@  ?  ok(1, $msg) : ok (0, join "\n", $@,  $msg );
 207 }
 208
 209 42;