lib/Math/BigInt/t/bigintc.t + VMS + perl@16925
[p5sagit/p5-mst-13.2.git] / lib / Math / Complex.pm
index daf190e..400366c 100644 (file)
@@ -7,20 +7,28 @@
 
 package Math::Complex;
 
-$VERSION = "1.30";
-
 our($VERSION, @ISA, @EXPORT, %EXPORT_TAGS, $Inf);
 
+$VERSION = 1.34;
+
 BEGIN {
-    eval { require POSIX; import POSIX 'HUGE_VAL' };
-    if (exists &HUGE_VAL) {
-       $Inf = sprintf "%g", &HUGE_VAL;
-    } else {   
-       my $e = $!;
-       $Inf = CORE::exp(CORE::exp(30));
-       $! = $e; # Clear ERANGE.
+    unless ($^O eq 'unicosmk') {
+        my $e = $!;
+       # We do want an arithmetic overflow, Inf INF inf Infinity:.
+        undef $Inf unless eval <<'EOE' and $Inf =~ /^inf(?:inity)?$/i;
+         local $SIG{FPE} = sub {die};
+         my $t = CORE::exp 30;
+         $Inf = CORE::exp $t;
+EOE
+       if (!defined $Inf) {            # Try a different method
+         undef $Inf unless eval <<'EOE' and $Inf =~ /^inf(?:inity)?$/i;
+           local $SIG{FPE} = sub {die};
+           my $t = 1;
+           $Inf = $t + "1e99999999999999999999999999999999";
+EOE
+       }
+        $! = $e; # Clear ERANGE.
     }
-    undef $Inf unless $Inf =~ /^inf(?:inity)?$/i; # Inf INF inf Infinity
     $Inf = "Inf" if !defined $Inf || !($Inf > 0); # Desperation.
 }
 
@@ -29,6 +37,9 @@ use strict;
 my $i;
 my %LOGN;
 
+# Regular expression for floating point numbers.
+my $gre = qr'\s*([\+\-]?(?:(?:(?:\d+(?:_\d+)*(?:\.\d*(?:_\d+)*)?|\.\d+(?:_\d+)*)(?:[eE][\+\-]?\d+(?:_\d+)*)?)))';
+
 require Exporter;
 
 @ISA = qw(Exporter);
@@ -100,6 +111,26 @@ sub _cannot_make {
     die "@{[(caller(1))[3]]}: Cannot take $_[0] of $_[1].\n";
 }
 
+sub _remake {
+    my $arg = shift;
+    my ($made, $p, $q);
+
+    if ($arg =~ /^(?:$gre)?$gre\s*i\s*$/) {
+       ($p, $q) = ($1 || 0, $2);
+       $made = 'cart';
+    } elsif ($arg =~ /^\s*\[\s*$gre\s*(?:,\s*$gre\s*)?\]\s*$/) {
+       ($p, $q) = ($1, $2 || 0);
+       $made = 'exp';
+    }
+
+    if ($made) {
+       $p =~ s/^\+//;
+       $q =~ s/^\+//;
+    }
+
+    return ($made, $p, $q);
+}
+
 #
 # ->make
 #
@@ -108,6 +139,16 @@ sub _cannot_make {
 sub make {
        my $self = bless {}, shift;
        my ($re, $im) = @_;
+       if (@_ == 1) {
+           my ($remade, $p, $q) = _remake($re);
+           if ($remade) {
+               if ($remade eq 'cart') {
+                   ($re, $im) = ($p, $q);
+               } else {
+                   return (ref $self)->emake($p, $q);
+               }
+           }
+       }
        my $rre = ref $re;
        if ( $rre ) {
            if ( $rre eq ref $self ) {
@@ -124,6 +165,9 @@ sub make {
                _cannot_make("imaginary part", $rim);
            }
        }
+       _cannot_make("real part",      $re) unless $re =~ /^$gre$/;
+       $im ||= 0;
+       _cannot_make("imaginary part", $im) unless $im =~ /^$gre$/;
        $self->{'cartesian'} = [ $re, $im ];
        $self->{c_dirty} = 0;
        $self->{p_dirty} = 1;
@@ -139,6 +183,16 @@ sub make {
 sub emake {
        my $self = bless {}, shift;
        my ($rho, $theta) = @_;
+       if (@_ == 1) {
+           my ($remade, $p, $q) = _remake($rho);
+           if ($remade) {
+               if ($remade eq 'exp') {
+                   ($rho, $theta) = ($p, $q);
+               } else {
+                   return (ref $self)->make($p, $q);
+               }
+           }
+       }
        my $rrh = ref $rho;
        if ( $rrh ) {
            if ( $rrh eq ref $self ) {
@@ -159,6 +213,9 @@ sub emake {
            $rho   = -$rho;
            $theta = ($theta <= 0) ? $theta + pi() : $theta - pi();
        }
+       _cannot_make("rho",   $rho)   unless $rho   =~ /^$gre$/;
+       $theta ||= 0;
+       _cannot_make("theta", $theta) unless $theta =~ /^$gre$/;
        $self->{'polar'} = [$rho, $theta];
        $self->{p_dirty} = 0;
        $self->{c_dirty} = 1;
@@ -175,8 +232,7 @@ sub new { &make }           # For backward compatibility only.
 # This avoids the burden of writing Math::Complex->make(re, im).
 #
 sub cplx {
-       my ($re, $im) = @_;
-       return __PACKAGE__->make($re, defined $im ? $im : 0);
+       return __PACKAGE__->make(@_);
 }
 
 #
@@ -186,8 +242,7 @@ sub cplx {
 # This avoids the burden of writing Math::Complex->emake(rho, theta).
 #
 sub cplxe {
-       my ($rho, $theta) = @_;
-       return __PACKAGE__->emake($rho, defined $theta ? $theta : 0);
+       return __PACKAGE__->emake(@_);
 }
 
 #
@@ -666,7 +721,7 @@ sub Re {
 #
 sub Im {
        my ($z, $Im) = @_;
-       return $z unless ref $z;
+       return 0 unless ref $z;
        if (defined $Im) {
            $z->{'cartesian'} = [ ${$z->cartesian}[0], $Im ];
            $z->{c_dirty} = 0;
@@ -879,7 +934,8 @@ sub acos {
        my $z = $_[0];
        return CORE::atan2(CORE::sqrt(1-$z*$z), $z)
            if (! ref $z) && CORE::abs($z) <= 1;
-       my ($x, $y) = ref $z ? @{$z->cartesian} : ($z, 0);
+       $z = cplx($z, 0) unless ref $z;
+       my ($x, $y) = @{$z->cartesian};
        return 0 if $x == 1 && $y == 0;
        my $t1 = CORE::sqrt(($x+1)*($x+1) + $y*$y);
        my $t2 = CORE::sqrt(($x-1)*($x-1) + $y*$y);
@@ -891,7 +947,7 @@ sub acos {
        my $u = CORE::atan2(CORE::sqrt(1-$beta*$beta), $beta);
        my $v = CORE::log($alpha + CORE::sqrt($alpha*$alpha-1));
        $v = -$v if $y > 0 || ($y == 0 && $x < -1);
-       return __PACKAGE__->make($u, $v);
+       return (ref $z)->make($u, $v);
 }
 
 #
@@ -903,7 +959,8 @@ sub asin {
        my $z = $_[0];
        return CORE::atan2($z, CORE::sqrt(1-$z*$z))
            if (! ref $z) && CORE::abs($z) <= 1;
-       my ($x, $y) = ref $z ? @{$z->cartesian} : ($z, 0);
+       $z = cplx($z, 0) unless ref $z;
+       my ($x, $y) = @{$z->cartesian};
        return 0 if $x == 0 && $y == 0;
        my $t1 = CORE::sqrt(($x+1)*($x+1) + $y*$y);
        my $t2 = CORE::sqrt(($x-1)*($x-1) + $y*$y);
@@ -915,7 +972,7 @@ sub asin {
        my $u =  CORE::atan2($beta, CORE::sqrt(1-$beta*$beta));
        my $v = -CORE::log($alpha + CORE::sqrt($alpha*$alpha-1));
        $v = -$v if $y > 0 || ($y == 0 && $x < -1);
-       return __PACKAGE__->make($u, $v);
+       return (ref $z)->make($u, $v);
 }
 
 #
@@ -1106,11 +1163,13 @@ sub acosh {
                if CORE::abs($re) < 1;
        }
        my $t = &sqrt($z * $z - 1) + $z;
-       # Try MacLaurin if looking bad.
+       # Try Taylor if looking bad (this usually means that
+       # $z was large negative, therefore the sqrt is really
+       # close to abs(z), summing that with z...)
        $t = 1/(2 * $z) - 1/(8 * $z**3) + 1/(16 * $z**5) - 5/(128 * $z**7)
            if $t == 0;
        my $u = &log($t);
-       $u->Im(-$u->Im) if $im == 0;
+       $u->Im(-$u->Im) if $re < 0 && $im == 0;
        return $re < 0 ? -$u : $u;
 }
 
@@ -1126,7 +1185,9 @@ sub asinh {
            return CORE::log($t) if $t;
        }
        my $t = &sqrt($z * $z + 1) + $z;
-       # Try MacLaurin if looking bad.
+       # Try Taylor if looking bad (this usually means that
+       # $z was large negative, therefore the sqrt is really
+       # close to abs(z), summing that with z...)
        $t = 1/(2 * $z) - 1/(8 * $z**3) + 1/(16 * $z**5) - 5/(128 * $z**7)
            if $t == 0;
        return &log($t);
@@ -1249,23 +1310,15 @@ sub display_format {
                my %obj = %{$self->{display_format}};
                @display_format{keys %obj} = values %obj;
            }
-           if (@_ == 1) {
-               $display_format{style} = shift;
-           } else {
-               my %new = @_;
-               @display_format{keys %new} = values %new;
-           }
-       } else {                                # Called as a class method
-           if (@_ = 1) {
-               $display_format{style} = $self;
-           } else {
-               my %new = @_;
-               @display_format{keys %new} = values %new;
-           }
-           undef $self;
+       }
+       if (@_ == 1) {
+           $display_format{style} = shift;
+       } else {
+           my %new = @_;
+           @display_format{keys %new} = values %new;
        }
 
-       if (defined $self) {
+       if (ref $self) { # Called as an object method
            $self->{display_format} = { %display_format };
            return
                wantarray ?
@@ -1273,6 +1326,7 @@ sub display_format {
                    $self->{display_format}->{style};
        }
 
+        # Called as a class method
        %DISPLAY_FORMAT = %display_format;
        return
            wantarray ?
@@ -1329,15 +1383,16 @@ sub stringify_cartesian {
        }
 
        if ($y) {
-           if ($y == 1)     { $im = ""  }
-           elsif ($y == -1) { $im = "-" }
-           elsif ($y =~ /^(NaN[QS]?)$/i) {
+           if ($y =~ /^(NaN[QS]?)$/i) {
                $im = $y;
            } else {
                if ($y =~ /^-?$Inf$/oi) {
                    $im = $y;
                } else {
-                   $im = defined $format ? sprintf($format, $y) : $y;
+                   $im =
+                       defined $format ?
+                           sprintf($format, $y) :
+                           ($y == 1 ? "" : ($y == -1 ? "-" : $y));
                }
            }
            $im .= "i";
@@ -1391,11 +1446,11 @@ sub stringify_polar {
 
        $t -= int(CORE::abs($t) / pit2) * pit2;
 
-       if ($format{polar_pretty_print}) {
+       if ($format{polar_pretty_print} && $t) {
            my ($a, $b);
            for $a (2..9) {
                $b = $t * $a / pi;
-               if (int($b) == $b) {
+               if ($b =~ /^-?\d+$/) {
                    $b = $b < 0 ? "-" : "" if CORE::abs($b) == 1;
                    $theta = "${b}pi/$a";
                    last;
@@ -1416,6 +1471,8 @@ sub stringify_polar {
 1;
 __END__
 
+=pod
+
 =head1 NAME
 
 Math::Complex - complex numbers and associated mathematical functions
@@ -1551,7 +1608,7 @@ be called an extension, would it?).
 
 A I<new> operation possible on a complex number that is
 the identity for real numbers is called the I<conjugate>, and is noted
-with an horizontal bar above the number, or C<~z> here.
+with a horizontal bar above the number, or C<~z> here.
 
         z = a + bi
        ~z = a - bi
@@ -1650,7 +1707,7 @@ I<arg>, I<abs>, I<log>, I<csc>, I<cot>, I<acsc>, I<acot>, I<csch>,
 I<coth>, I<acosech>, I<acotanh>, have aliases I<rho>, I<theta>, I<ln>,
 I<cosec>, I<cotan>, I<acosec>, I<acotan>, I<cosech>, I<cotanh>,
 I<acosech>, I<acotanh>, respectively.  C<Re>, C<Im>, C<arg>, C<abs>,
-C<rho>, and C<theta> can be used also also mutators.  The C<cbrt>
+C<rho>, and C<theta> can be used also as mutators.  The C<cbrt>
 returns only one of the solutions: if you want all three, use the
 C<root> function.
 
@@ -1703,13 +1760,25 @@ but that will be silently converted into C<[3,-3pi/4]>, since the
 modulus must be non-negative (it represents the distance to the origin
 in the complex plane).
 
-It is also possible to have a complex number as either argument of
-either the C<make> or C<emake>: the appropriate component of
+It is also possible to have a complex number as either argument of the
+C<make>, C<emake>, C<cplx>, and C<cplxe>: the appropriate component of
 the argument will be used.
 
        $z1 = cplx(-2,  1);
        $z2 = cplx($z1, 4);
 
+The C<new>, C<make>, C<emake>, C<cplx>, and C<cplxe> will also
+understand a single (string) argument of the forms
+
+       2-3i
+       -3i
+       [2,3]
+       [2]
+
+in which case the appropriate cartesian and exponential components
+will be parsed from the string and used to create new complex numbers.
+The imaginary component and the theta, respectively, will default to zero.
+
 =head1 STRINGIFICATION
 
 When printed, a complex number is usually shown under its cartesian
@@ -1747,29 +1816,33 @@ C<display_format> object method can now be called using
 a parameter hash instead of just a one parameter.
 
 The old display format style, which can have values C<"cartesian"> or
-C<"polar">, can be changed using the C<"style"> parameter.  (The one
-parameter calling convention also still works.)
+C<"polar">, can be changed using the C<"style"> parameter.
+
+       $j->display_format(style => "polar");
+
+The one parameter calling convention also still works.
+
+       $j->display_format("polar");
 
 There are two new display parameters.
 
-The first one is C<"format">, which is a sprintf()-style format
-string to be used for both parts of the complex number(s).  The
-default is C<undef>, which corresponds usually (this is somewhat
-system-dependent) to C<"%.15g">.  You can revert to the default by
-setting the format string to C<undef>.
+The first one is C<"format">, which is a sprintf()-style format string
+to be used for both numeric parts of the complex number(s).  The is
+somewhat system-dependent but most often it corresponds to C<"%.15g">.
+You can revert to the default by setting the C<format> to C<undef>.
 
        # the $j from the above example
 
        $j->display_format('format' => '%.5f');
        print "j = $j\n";               # Prints "j = -0.50000+0.86603i"
-       $j->display_format('format' => '%.6f');
+       $j->display_format('format' => undef);
        print "j = $j\n";               # Prints "j = -0.5+0.86603i"
 
 Notice that this affects also the return values of the
 C<display_format> methods: in list context the whole parameter hash
-will be returned, as opposed to only the style parameter value.  If
-you want to know the whole truth for a complex number, you must call
-both the class method and the object method:
+will be returned, as opposed to only the style parameter value.
+This is a potential incompatibility with earlier versions if you
+have been calling the C<display_format> method in list context.
 
 The second new display parameter is C<"polar_pretty_print">, which can
 be set to true or false, the default being true.  See the previous
@@ -1822,7 +1895,7 @@ or
        Died at...
 
 For the C<csc>, C<cot>, C<asec>, C<acsc>, C<acot>, C<csch>, C<coth>,
-C<asech>, C<acsch>, the argument cannot be C<0> (zero).  For the the
+C<asech>, C<acsch>, the argument cannot be C<0> (zero).  For the
 logarithmic functions and the C<atanh>, C<acoth>, the argument cannot
 be C<1> (one).  For the C<atanh>, C<acoth>, the argument cannot be
 C<-1> (minus one).  For the C<atan>, C<acot>, the argument cannot be
@@ -1833,8 +1906,7 @@ is any integer.
 
 Note that because we are operating on approximations of real numbers,
 these errors can happen when merely `too close' to the singularities
-listed above.  For example C<tan(2*atan2(1,1)+1e-15)> will die of
-division by zero.
+listed above.
 
 =head1 ERRORS DUE TO INDIGESTIBLE ARGUMENTS
 
@@ -1864,10 +1936,10 @@ Whatever it is, it does not manifest itself anywhere else where Perl runs.
 
 =head1 AUTHORS
 
-Raphael Manfredi <F<Raphael_Manfredi@pobox.com>> and
-Jarkko Hietaniemi <F<jhi@iki.fi>>.
+Daniel S. Lewart <F<d-lewart@uiuc.edu>>
 
-Extensive patches by Daniel S. Lewart <F<d-lewart@uiuc.edu>>.
+Original authors Raphael Manfredi <F<Raphael_Manfredi@pobox.com>> and
+Jarkko Hietaniemi <F<jhi@iki.fi>>
 
 =cut