blob: 5ac6b186c6a8abbed15f1b445b0fdbbf60da329d [file] [log] [blame]
Guido van Rossume8769491992-08-13 12:14:11 +00001# Complex numbers
Guido van Rossum81a12bc1994-10-08 18:56:41 +00002# ---------------
3
4# This module represents complex numbers as instances of the class Complex.
5# A Complex instance z has two data attribues, z.re (the real part) and z.im
6# (the imaginary part). In fact, z.re and z.im can have any value -- all
7# arithmetic operators work regardless of the type of z.re and z.im (as long
8# as they support numerical operations).
9#
10# The following functions exist (Complex is actually a class):
11# Complex([re [,im]) -> creates a complex number from a real and an imaginary part
12# IsComplex(z) -> true iff z is a complex number (== has .re and .im attributes)
13# ToComplex(z) -> a complex number equal to z; z itself if IsComplex(z) is true
14# if z is a tuple(re, im) it will also be converted
15# PolarToComplex([r [,phi [,fullcircle]]]) ->
16# the complex number z for which r == z.radius() and phi == z.angle(fullcircle)
17# (r and phi default to 0)
18#
19# Complex numbers have the following methods:
20# z.abs() -> absolute value of z
21# z.radius() == z.abs()
22# z.angle([fullcircle]) -> angle from positive X axis; fullcircle gives units
23# z.phi([fullcircle]) == z.angle(fullcircle)
24#
25# These standard functions and unary operators accept complex arguments:
26# abs(z)
27# -z
28# +z
29# not z
30# repr(z) == `z`
31# str(z)
32# hash(z) -> a combination of hash(z.re) and hash(z.im) such that if z.im is zero
33# the result equals hash(z.re)
34# Note that hex(z) and oct(z) are not defined.
35#
36# These conversions accept complex arguments only if their imaginary part is zero:
37# int(z)
38# long(z)
39# float(z)
40#
41# The following operators accept two complex numbers, or one complex number
42# and one real number (int, long or float):
43# z1 + z2
44# z1 - z2
45# z1 * z2
46# z1 / z2
47# pow(z1, z2)
48# cmp(z1, z2)
49# Note that z1 % z2 and divmod(z1, z2) are not defined,
50# nor are shift and mask operations.
51#
52# The standard module math does not support complex numbers.
53# (I suppose it would be easy to implement a cmath module.)
54#
55# Idea:
56# add a class Polar(r, phi) and mixed-mode arithmetic which
57# chooses the most appropriate type for the result:
58# Complex for +,-,cmp
59# Polar for *,/,pow
Guido van Rossume8769491992-08-13 12:14:11 +000060
61
Guido van Rossum81a12bc1994-10-08 18:56:41 +000062import types, math
Guido van Rossume8769491992-08-13 12:14:11 +000063
Guido van Rossum81a12bc1994-10-08 18:56:41 +000064twopi = math.pi*2.0
65halfpi = math.pi/2.0
Guido van Rossume8769491992-08-13 12:14:11 +000066
Guido van Rossum81a12bc1994-10-08 18:56:41 +000067def IsComplex(obj):
68 return hasattr(obj, 're') and hasattr(obj, 'im')
Guido van Rossume8769491992-08-13 12:14:11 +000069
Guido van Rossum81a12bc1994-10-08 18:56:41 +000070def ToComplex(obj):
71 if IsComplex(obj):
72 return obj
73 elif type(obj) == types.TupleType:
74 return apply(Complex, obj)
75 else:
76 return Complex(obj)
Guido van Rossum7565b931993-12-17 14:23:52 +000077
Guido van Rossum81a12bc1994-10-08 18:56:41 +000078def PolarToComplex(r = 0, phi = 0, fullcircle = twopi):
79 phi = phi * (twopi / fullcircle)
80 return Complex(math.cos(phi)*r, math.sin(phi)*r)
81
82def Re(obj):
83 if IsComplex(obj):
84 return obj.re
85 else:
86 return obj
87
88def Im(obj):
89 if IsComplex(obj):
90 return obj.im
91 else:
92 return obj
93
94class Complex:
95
96 def __init__(self, re=0, im=0):
97 if IsComplex(re):
98 im = i + Complex(0, re.im)
99 re = re.re
100 if IsComplex(im):
101 re = re - im.im
102 im = im.re
103 self.__dict__['re'] = re
104 self.__dict__['im'] = im
105
106 def __setattr__(self, name, value):
107 raise TypeError, 'Complex numbers are immutable'
108
109 def __hash__(self):
110 if not self.im: return hash(self.re)
111 mod = sys.maxint + 1L
112 return int((hash(self.re) + 2L*hash(self.im) + mod) % (2L*mod) - mod)
Guido van Rossume8769491992-08-13 12:14:11 +0000113
114 def __repr__(self):
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000115 if not self.im:
116 return 'Complex(%s)' % `self.re`
117 else:
118 return 'Complex(%s, %s)' % (`self.re`, `self.im`)
Guido van Rossume8769491992-08-13 12:14:11 +0000119
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000120 def __str__(self):
121 if not self.im:
122 return `self.re`
123 else:
124 return 'Complex(%s, %s)' % (`self.re`, `self.im`)
125
126 def __neg__(self):
127 return Complex(-self.re, -self.im)
128
129 def __pos__(self):
130 return self
131
132 def __abs__(self):
133 # XXX could be done differently to avoid overflow!
134 return math.sqrt(self.re*self.re + self.im*self.im)
135
136 def __int__(self):
137 if self.im:
138 raise ValueError, "can't convert Complex with nonzero im to int"
139 return int(self.re)
140
141 def __long__(self):
142 if self.im:
143 raise ValueError, "can't convert Complex with nonzero im to long"
144 return long(self.re)
Guido van Rossume8769491992-08-13 12:14:11 +0000145
146 def __float__(self):
147 if self.im:
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000148 raise ValueError, "can't convert Complex with nonzero im to float"
Guido van Rossume8769491992-08-13 12:14:11 +0000149 return float(self.re)
150
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000151 def __cmp__(self, other):
152 other = ToComplex(other)
153 return cmp((self.re, self.im), (other.re, other.im))
Guido van Rossume8769491992-08-13 12:14:11 +0000154
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000155 def __rcmp__(self, other):
156 other = ToComplex(other)
157 return cmp(other, self)
158
159 def __nonzero__(self):
160 return not (self.re == self.im == 0)
Guido van Rossume8769491992-08-13 12:14:11 +0000161
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000162 abs = radius = __abs__
Guido van Rossume8769491992-08-13 12:14:11 +0000163
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000164 def angle(self, fullcircle = twopi):
165 return (fullcircle/twopi) * ((halfpi - math.atan2(self.re, self.im)) % twopi)
Guido van Rossume8769491992-08-13 12:14:11 +0000166
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000167 phi = angle
Guido van Rossume8769491992-08-13 12:14:11 +0000168
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000169 def __add__(self, other):
170 other = ToComplex(other)
171 return Complex(self.re + other.re, self.im + other.im)
Guido van Rossume8769491992-08-13 12:14:11 +0000172
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000173 __radd__ = __add__
Guido van Rossume8769491992-08-13 12:14:11 +0000174
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000175 def __sub__(self, other):
176 other = ToComplex(other)
177 return Complex(self.re - other.re, self.im - other.im)
178
179 def __rsub__(self, other):
180 other = ToComplex(other)
181 return other - self
182
183 def __mul__(self, other):
184 other = ToComplex(other)
185 return Complex(self.re*other.re - self.im*other.im,
186 self.re*other.im + self.im*other.re)
187
188 __rmul__ = __mul__
189
190 def __div__(self, other):
191 other = ToComplex(other)
192 d = float(other.re*other.re + other.im*other.im)
193 if not d: raise ZeroDivisionError, 'Complex division'
194 return Complex((self.re*other.re + self.im*other.im) / d,
195 (self.im*other.re - self.re*other.im) / d)
196
197 def __rdiv__(self, other):
198 other = ToComplex(other)
199 return other / self
200
201 def __pow__(self, n, z=None):
202 if z is not None:
203 raise TypeError, 'Complex does not support ternary pow()'
204 if IsComplex(n):
205 if n.im: raise TypeError, 'Complex to the Complex power'
206 n = n.re
207 r = pow(self.abs(), n)
208 phi = n*self.angle()
209 return Complex(math.cos(phi)*r, math.sin(phi)*r)
210
211 def __rpow__(self, base):
212 base = ToComplex(base)
213 return pow(base, self)
214
215
216def checkop(expr, a, b, value, fuzz = 1e-6):
217 import sys
218 print ' ', a, 'and', b,
219 try:
220 result = eval(expr)
221 except:
222 result = sys.exc_type
223 print '->', result
224 if (type(result) == type('') or type(value) == type('')):
225 ok = result == value
226 else:
227 ok = abs(result - value) <= fuzz
228 if not ok:
229 print '!!\t!!\t!! should be', value, 'diff', abs(result - value)
Guido van Rossume8769491992-08-13 12:14:11 +0000230
231
232def test():
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000233 testsuite = {
234 'a+b': [
235 (1, 10, 11),
236 (1, Complex(0,10), Complex(1,10)),
237 (Complex(0,10), 1, Complex(1,10)),
238 (Complex(0,10), Complex(1), Complex(1,10)),
239 (Complex(1), Complex(0,10), Complex(1,10)),
240 ],
241 'a-b': [
242 (1, 10, -9),
243 (1, Complex(0,10), Complex(1,-10)),
244 (Complex(0,10), 1, Complex(-1,10)),
245 (Complex(0,10), Complex(1), Complex(-1,10)),
246 (Complex(1), Complex(0,10), Complex(1,-10)),
247 ],
248 'a*b': [
249 (1, 10, 10),
250 (1, Complex(0,10), Complex(0, 10)),
251 (Complex(0,10), 1, Complex(0,10)),
252 (Complex(0,10), Complex(1), Complex(0,10)),
253 (Complex(1), Complex(0,10), Complex(0,10)),
254 ],
255 'a/b': [
256 (1., 10, 0.1),
257 (1, Complex(0,10), Complex(0, -0.1)),
258 (Complex(0, 10), 1, Complex(0, 10)),
259 (Complex(0, 10), Complex(1), Complex(0, 10)),
260 (Complex(1), Complex(0,10), Complex(0, -0.1)),
261 ],
262 'pow(a,b)': [
263 (1, 10, 1),
264 (1, Complex(0,10), 'TypeError'),
265 (Complex(0,10), 1, Complex(0,10)),
266 (Complex(0,10), Complex(1), Complex(0,10)),
267 (Complex(1), Complex(0,10), 'TypeError'),
268 (2, Complex(4,0), 16),
269 ],
270 'cmp(a,b)': [
271 (1, 10, -1),
272 (1, Complex(0,10), 1),
273 (Complex(0,10), 1, -1),
274 (Complex(0,10), Complex(1), -1),
275 (Complex(1), Complex(0,10), 1),
276 ],
277 }
278 exprs = testsuite.keys()
279 exprs.sort()
280 for expr in exprs:
281 print expr + ':'
282 t = (expr,)
283 for item in testsuite[expr]:
284 apply(checkop, t+item)
285
Guido van Rossume8769491992-08-13 12:14:11 +0000286
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000287if __name__ == '__main__':
288 test()