lib/Math/Complex.pm

   1 package Math::Complex;
   2
   3 require Exporter;
   4
   5 @ISA = ('Exporter');
   6
   7 # just to make use happy
   8
   9 use overload
  10     '+'   => sub  { my($x1,$y1,$x2,$y2) = (@{$_[0]},@{$_[1]});
  11                       bless [ $x1+$x2, $y1+$y2];
  12              },
  13
  14     '-'   => sub  { my($x1,$y1,$x2,$y2) = (@{$_[0]},@{$_[1]});
  15                       bless [ $x1-$x2, $y1-$y2];
  16              },
  17
  18     '*'   => sub  { my($x1,$y1,$x2,$y2) = (@{$_[0]},@{$_[1]});
  19                     bless [ $x1*$x2-$y1*$y2,$x1*$y2+$x2*$y1];
  20              },
  21
  22     '/'   => sub  { my($x1,$y1,$x2,$y2) = (@{$_[0]},@{$_[1]});
  23                     my $q = $x2*$x2+$y2*$y2;
  24                     bless [($x1*$x2+$y1*$y2)/$q, ($y1*$x2-$y2*$x1)/$q];
  25              },
  26
  27     'neg' => sub  { my($x,$y) = @{$_[0]}; bless [ -$x, -$y];
  28              },
  29
  30     '~'   => sub  { my($x,$y) = @{$_[0]}; bless [ $x, -$y];
  31              },
  32
  33     'abs'   => sub  { my($x,$y) = @{$_[0]}; sqrt $x*$x+$y*$y;
  34              },
  35
  36     'cos' => sub { my($x,$y) = @{$_[0]};
  37                    my ($ab,$c,$s) = (exp $y, cos $x, sin $x);
  38                    my $abr = 1/(2*$ab); $ab /= 2;
  39                    bless [ ($abr+$ab)*$c, ($abr-$ab)*$s];
  40              },
  41
  42     'sin' => sub { my($x,$y) = @{$_[0]};
  43                    my ($ab,$c,$s) = (exp $y, cos $x, sin $x);
  44                    my $abr = 1/(2*$ab); $ab /= 2;
  45                    bless [ (-$abr-$ab)*$s, ($abr-$ab)*$c];
  46              },
  47
  48     'exp' => sub { my($x,$y) = @{$_[0]};
  49                    my ($ab,$c,$s) = (exp $x, cos $y, sin $y);
  50                    bless [ $ab*$c, $ab*$s ];
  51               },
  52
  53     'sqrt' => sub {
  54         my($zr,$zi) = @{$_[0]};
  55         my ($x, $y, $r, $w);
  56         my $c = new Math::Complex (0,0);
  57         if (($zr == 0) && ($zi == 0)) {
  58             # nothing, $c already set
  59         }
  60         else {
  61           $x = abs($zr);
  62           $y = abs($zi);
  63           if ($x >= $y) {
  64               $r = $y/$x;
  65               $w = sqrt(0.5 * $x * (1.0+sqrt(1.0+$r*$r)));
  66           }
  67           else {
  68               $r = $x/$y;
  69               $w = sqrt(0.5 * ($x + $y*sqrt(1.0+$r*$r)));
  70           }
  71           if ( $zr >= 0) {
  72               @$c = ($w, $zi/(2 * $w) );
  73           }
  74           else {
  75               $c->[1] = ($zi >= 0) ? $w : -$w;
  76               $c->[0] = $zi/(2.0* $c->[1]);
  77           }
  78         }
  79         return $c;
  80       },
  81
  82     qw("" stringify)
  83 ;
  84
  85 sub new {
  86     my $class = shift;
  87     my @C = @_;
  88     bless \@C, $class;
  89 }
  90
  91 sub Re {
  92     my($x,$y) = @{$_[0]};
  93     $x;
  94 }
  95
  96 sub Im {
  97     my($x,$y) = @{$_[0]};
  98     $y;
  99 }
 100
 101 sub arg {
 102     my($x,$y) = @{$_[0]};
 103     atan2($y,$x);
 104 }
 105
 106 sub stringify {
 107     my($x,$y) = @{$_[0]};
 108     my($re,$im);
 109
 110     $re = $x if ($x);
 111     if ($y == 1) {$im = 'i';}
 112     elsif ($y == -1){$im = '-i';}
 113     elsif ($y) {$im = $y . 'i'; }
 114
 115     local $_ = $re.'+'.$im;
 116     s/\+-/-/;
 117     s/^\+//;
 118     s/[\+-]$//;
 119     $_ = 0 if ($_ eq '');
 120     return $_;
 121 }
 122
 123 1;
 124 __END__
 125
 126 =head1 NAME
 127
 128 Math::Complex - complex numbers package
 129
 130 =head1 SYNOPSIS
 131
 132   use Math::Complex;
 133   $i = new Math::Complex;
 134
 135 =head1 DESCRIPTION
 136
 137 Complex numbers declared as
 138
 139     $i = Math::Complex->new(1,1);
 140
 141 can be manipulated with overloaded math operators. The operators
 142
 143   + - * / neg ~ abs cos sin exp sqrt
 144
 145 are supported as well as
 146
 147   "" (stringify)
 148
 149 The methods
 150
 151   Re Im arg
 152
 153 are also provided.
 154
 155 =head1 BUGS
 156
 157 sqrt() should return two roots, but only returns one.
 158
 159 =head1 AUTHORS
 160
 161 Dave Nadler, Tom Christiansen, Tim Bunce, Larry Wall.
 162
 163 =cut