Commit | Line | Data |
a0d0e21e |
1 | # |
2 | # Perl5 Package for complex numbers |
3 | # |
4 | # 1994 by David Nadler |
5 | # Coding know-how provided by Tom Christiansen, Tim Bunce, and Larry Wall |
6 | # sqrt() added by Tom Christiansen; beware should have two roots, |
7 | # but only returns one. (use wantarray?) |
8 | # |
9 | # |
10 | # The functions "Re", "Im", and "arg" are provided. |
11 | # "~" is used as the conjugation operator and "abs" is overloaded. |
12 | # |
13 | # Transcendental functions overloaded: so far only sin, cos, and exp. |
14 | # |
15 | |
16 | package Math::Complex; |
17 | |
18 | require Exporter; |
19 | |
20 | @ISA = ('Exporter'); |
21 | |
22 | # just to make use happy |
23 | |
24 | %OVERLOAD= ( |
25 | '+' => sub { my($x1,$y1,$x2,$y2) = (@{$_[0]},@{$_[1]}); |
26 | bless [ $x1+$x2, $y1+$y2]; |
27 | }, |
28 | |
29 | '-' => sub { my($x1,$y1,$x2,$y2) = (@{$_[0]},@{$_[1]}); |
30 | bless [ $x1-$x2, $y1-$y2]; |
31 | }, |
32 | |
33 | '*' => sub { my($x1,$y1,$x2,$y2) = (@{$_[0]},@{$_[1]}); |
34 | bless [ $x1*$x2-$y1*$y2,$x1*$y2+$x2*$y1]; |
35 | }, |
36 | |
37 | '/' => sub { my($x1,$y1,$x2,$y2) = (@{$_[0]},@{$_[1]}); |
38 | my $q = $x2*$x2+$y2*$y2; |
39 | bless [($x1*$x2+$y1*$y2)/$q, ($y1*$x2-$y2*$x1)/$q]; |
40 | }, |
41 | |
42 | 'neg' => sub { my($x,$y) = @{$_[0]}; bless [ -$x, -$y]; |
43 | }, |
44 | |
45 | '~' => sub { my($x,$y) = @{$_[0]}; bless [ $x, -$y]; |
46 | }, |
47 | |
48 | 'abs' => sub { my($x,$y) = @{$_[0]}; sqrt $x*$x+$y*$y; |
49 | }, |
50 | |
51 | 'cos' => sub { my($x,$y) = @{$_[0]}; |
52 | my ($ab,$c,$s) = (exp $y, cos $x, sin $x); |
53 | my $abr = 1/(2*$ab); $ab /= 2; |
54 | bless [ ($abr+$ab)*$c, ($abr-$ab)*$s]; |
55 | }, |
56 | |
57 | 'sin' => sub { my($x,$y) = @{$_[0]}; |
58 | my ($ab,$c,$s) = (exp $y, cos $x, sin $x); |
59 | my $abr = 1/(2*$ab); $ab /= 2; |
60 | bless [ (-$abr-$ab)*$s, ($abr-$ab)*$c]; |
61 | }, |
62 | |
63 | 'exp' => sub { my($x,$y) = @{$_[0]}; |
64 | my ($ab,$c,$s) = (exp $x, cos $y, sin $y); |
65 | bless [ $ab*$c, $ab*$s ]; |
66 | }, |
67 | |
68 | 'sqrt' => sub { |
69 | my($zr,$zi) = @{$_[0]}; |
70 | my ($x, $y, $r, $w); |
71 | my $c = new Math::Complex (0,0); |
72 | if (($zr == 0) && ($zi == 0)) { |
73 | # nothing, $c already set |
74 | } |
75 | else { |
76 | $x = abs($zr); |
77 | $y = abs($zi); |
78 | if ($x >= $y) { |
79 | $r = $y/$x; |
80 | $w = sqrt($x) * sqrt(0.5*(1.0+sqrt(1.0+$r*$r))); |
81 | } |
82 | else { |
83 | $r = $x/$y; |
84 | $w = sqrt($y) * sqrt($y) * sqrt(0.5*($r+sqrt(1.0+$r*$r))); |
85 | } |
86 | if ( $zr >= 0) { |
87 | @$c = ($w, $zi/(2 * $w) ); |
88 | } |
89 | else { |
90 | $c->[1] = ($zi >= 0) ? $w : -$w; |
91 | $c->[0] = $zi/(2.0* $c->[1]); |
92 | } |
93 | } |
94 | return $c; |
95 | }, |
96 | |
97 | qw("" stringify) |
98 | ); |
99 | |
100 | sub new { |
101 | shift; |
102 | my @C = @_; |
103 | bless \@C; |
104 | } |
105 | |
106 | sub Re { |
107 | my($x,$y) = @{$_[0]}; |
108 | $x; |
109 | } |
110 | |
111 | sub Im { |
112 | my($x,$y) = @{$_[0]}; |
113 | $y; |
114 | } |
115 | |
116 | sub arg { |
117 | my($x,$y) = @{$_[0]}; |
118 | atan2($y,$x); |
119 | } |
120 | |
121 | sub stringify { |
122 | my($x,$y) = @{$_[0]}; |
123 | my($re,$im); |
124 | |
125 | $re = $x if ($x); |
126 | if ($y == 1) {$im = 'i';} |
127 | elsif ($y == -1){$im = '-i';} |
128 | elsif ($y) {$im = "${y}i"; } |
129 | |
130 | local $_ = $re.'+'.$im; |
131 | s/\+-/-/; |
132 | s/^\+//; |
133 | s/[\+-]$//; |
134 | $_ = 0 if ($_ eq ''); |
135 | return $_; |
136 | } |