blob: 06f2212afa611e9c183d98b5885ca192aed62d59 [file] [log] [blame]
Guido van Rossume8769491992-08-13 12:14:11 +00001# Complex numbers
Guido van Rossum81a12bc1994-10-08 18:56:41 +00002# ---------------
3
Guido van Rossum72ba6161996-07-30 19:02:01 +00004# [Now that Python has a complex data type built-in, this is not very
5# useful, but it's still a nice example class]
6
Guido van Rossum81a12bc1994-10-08 18:56:41 +00007# This module represents complex numbers as instances of the class Complex.
8# A Complex instance z has two data attribues, z.re (the real part) and z.im
9# (the imaginary part). In fact, z.re and z.im can have any value -- all
10# arithmetic operators work regardless of the type of z.re and z.im (as long
11# as they support numerical operations).
12#
13# The following functions exist (Complex is actually a class):
14# Complex([re [,im]) -> creates a complex number from a real and an imaginary part
15# IsComplex(z) -> true iff z is a complex number (== has .re and .im attributes)
16# ToComplex(z) -> a complex number equal to z; z itself if IsComplex(z) is true
17# if z is a tuple(re, im) it will also be converted
18# PolarToComplex([r [,phi [,fullcircle]]]) ->
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +000019# the complex number z for which r == z.radius() and phi == z.angle(fullcircle)
20# (r and phi default to 0)
Guido van Rossum72ba6161996-07-30 19:02:01 +000021# exp(z) -> returns the complex exponential of z. Equivalent to pow(math.e,z).
Guido van Rossum81a12bc1994-10-08 18:56:41 +000022#
23# Complex numbers have the following methods:
24# z.abs() -> absolute value of z
25# z.radius() == z.abs()
26# z.angle([fullcircle]) -> angle from positive X axis; fullcircle gives units
27# z.phi([fullcircle]) == z.angle(fullcircle)
28#
29# These standard functions and unary operators accept complex arguments:
30# abs(z)
31# -z
32# +z
33# not z
34# repr(z) == `z`
35# str(z)
36# hash(z) -> a combination of hash(z.re) and hash(z.im) such that if z.im is zero
37# the result equals hash(z.re)
38# Note that hex(z) and oct(z) are not defined.
39#
40# These conversions accept complex arguments only if their imaginary part is zero:
41# int(z)
Guido van Rossum81a12bc1994-10-08 18:56:41 +000042# float(z)
43#
44# The following operators accept two complex numbers, or one complex number
45# and one real number (int, long or float):
46# z1 + z2
47# z1 - z2
48# z1 * z2
49# z1 / z2
50# pow(z1, z2)
51# cmp(z1, z2)
52# Note that z1 % z2 and divmod(z1, z2) are not defined,
53# nor are shift and mask operations.
54#
55# The standard module math does not support complex numbers.
Raymond Hettingerbdaad8c2005-04-09 14:55:07 +000056# The cmath modules should be used instead.
Guido van Rossum81a12bc1994-10-08 18:56:41 +000057#
58# Idea:
59# add a class Polar(r, phi) and mixed-mode arithmetic which
60# chooses the most appropriate type for the result:
61# Complex for +,-,cmp
62# Polar for *,/,pow
Guido van Rossume8769491992-08-13 12:14:11 +000063
Raymond Hettingerbdaad8c2005-04-09 14:55:07 +000064import math
Martin v. Löwis4a1e48c2005-04-09 10:51:19 +000065import sys
Guido van Rossume8769491992-08-13 12:14:11 +000066
Guido van Rossum81a12bc1994-10-08 18:56:41 +000067twopi = math.pi*2.0
68halfpi = math.pi/2.0
Guido van Rossume8769491992-08-13 12:14:11 +000069
Guido van Rossum81a12bc1994-10-08 18:56:41 +000070def IsComplex(obj):
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +000071 return hasattr(obj, 're') and hasattr(obj, 'im')
Guido van Rossume8769491992-08-13 12:14:11 +000072
Guido van Rossum81a12bc1994-10-08 18:56:41 +000073def ToComplex(obj):
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +000074 if IsComplex(obj):
75 return obj
Raymond Hettingerbdaad8c2005-04-09 14:55:07 +000076 elif isinstance(obj, tuple):
77 return Complex(*obj)
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +000078 else:
79 return Complex(obj)
Guido van Rossum7565b931993-12-17 14:23:52 +000080
Guido van Rossum81a12bc1994-10-08 18:56:41 +000081def PolarToComplex(r = 0, phi = 0, fullcircle = twopi):
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +000082 phi = phi * (twopi / fullcircle)
83 return Complex(math.cos(phi)*r, math.sin(phi)*r)
Guido van Rossum81a12bc1994-10-08 18:56:41 +000084
85def Re(obj):
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +000086 if IsComplex(obj):
87 return obj.re
Raymond Hettingerbdaad8c2005-04-09 14:55:07 +000088 return obj
Guido van Rossum81a12bc1994-10-08 18:56:41 +000089
90def Im(obj):
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +000091 if IsComplex(obj):
92 return obj.im
Raymond Hettingerbdaad8c2005-04-09 14:55:07 +000093 return 0
Guido van Rossum81a12bc1994-10-08 18:56:41 +000094
95class Complex:
96
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +000097 def __init__(self, re=0, im=0):
Martin v. Löwis4a1e48c2005-04-09 10:51:19 +000098 _re = 0
99 _im = 0
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000100 if IsComplex(re):
Martin v. Löwis4a1e48c2005-04-09 10:51:19 +0000101 _re = re.re
102 _im = re.im
103 else:
104 _re = re
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000105 if IsComplex(im):
Martin v. Löwis4a1e48c2005-04-09 10:51:19 +0000106 _re = _re - im.im
107 _im = _im + im.re
108 else:
109 _im = _im + im
110 # this class is immutable, so setting self.re directly is
111 # not possible.
112 self.__dict__['re'] = _re
113 self.__dict__['im'] = _im
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000114
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000115 def __setattr__(self, name, value):
Collin Winter6f2df4d2007-07-17 20:59:35 +0000116 raise TypeError('Complex numbers are immutable')
Guido van Rossume8769491992-08-13 12:14:11 +0000117
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000118 def __hash__(self):
Raymond Hettingerbdaad8c2005-04-09 14:55:07 +0000119 if not self.im:
120 return hash(self.re)
121 return hash((self.re, self.im))
Guido van Rossume8769491992-08-13 12:14:11 +0000122
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000123 def __repr__(self):
124 if not self.im:
Walter Dörwald70a6b492004-02-12 17:35:32 +0000125 return 'Complex(%r)' % (self.re,)
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000126 else:
Walter Dörwald70a6b492004-02-12 17:35:32 +0000127 return 'Complex(%r, %r)' % (self.re, self.im)
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000128
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000129 def __str__(self):
130 if not self.im:
Walter Dörwald70a6b492004-02-12 17:35:32 +0000131 return repr(self.re)
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000132 else:
Walter Dörwald70a6b492004-02-12 17:35:32 +0000133 return 'Complex(%r, %r)' % (self.re, self.im)
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000134
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000135 def __neg__(self):
136 return Complex(-self.re, -self.im)
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000137
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000138 def __pos__(self):
139 return self
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000140
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000141 def __abs__(self):
Raymond Hettingerbdaad8c2005-04-09 14:55:07 +0000142 return math.hypot(self.re, self.im)
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000143
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000144 def __int__(self):
145 if self.im:
Collin Winter6f2df4d2007-07-17 20:59:35 +0000146 raise ValueError("can't convert Complex with nonzero im to int")
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000147 return int(self.re)
Guido van Rossume8769491992-08-13 12:14:11 +0000148
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000149 def __float__(self):
150 if self.im:
Collin Winter6f2df4d2007-07-17 20:59:35 +0000151 raise ValueError("can't convert Complex with nonzero im to float")
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000152 return float(self.re)
Guido van Rossume8769491992-08-13 12:14:11 +0000153
Alexander Belopolsky74135d02010-07-03 22:36:06 +0000154 def __eq__(self, other):
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000155 other = ToComplex(other)
Alexander Belopolsky74135d02010-07-03 22:36:06 +0000156 return (self.re, self.im) == (other.re, other.im)
Guido van Rossume8769491992-08-13 12:14:11 +0000157
Jack Diederich4dafcc42006-11-28 19:15:13 +0000158 def __bool__(self):
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000159 return not (self.re == self.im == 0)
Guido van Rossume8769491992-08-13 12:14:11 +0000160
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000161 abs = radius = __abs__
Guido van Rossume8769491992-08-13 12:14:11 +0000162
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000163 def angle(self, fullcircle = twopi):
164 return (fullcircle/twopi) * ((halfpi - math.atan2(self.re, self.im)) % twopi)
Guido van Rossume8769491992-08-13 12:14:11 +0000165
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000166 phi = angle
Guido van Rossume8769491992-08-13 12:14:11 +0000167
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000168 def __add__(self, other):
169 other = ToComplex(other)
170 return Complex(self.re + other.re, self.im + other.im)
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000171
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000172 __radd__ = __add__
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000173
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000174 def __sub__(self, other):
175 other = ToComplex(other)
176 return Complex(self.re - other.re, self.im - other.im)
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000177
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000178 def __rsub__(self, other):
179 other = ToComplex(other)
180 return other - self
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000181
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000182 def __mul__(self, other):
183 other = ToComplex(other)
184 return Complex(self.re*other.re - self.im*other.im,
185 self.re*other.im + self.im*other.re)
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000186
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000187 __rmul__ = __mul__
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000188
Alexander Belopolsky74135d02010-07-03 22:36:06 +0000189 def __truediv__(self, other):
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000190 other = ToComplex(other)
191 d = float(other.re*other.re + other.im*other.im)
Collin Winter6f2df4d2007-07-17 20:59:35 +0000192 if not d: raise ZeroDivisionError('Complex division')
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000193 return Complex((self.re*other.re + self.im*other.im) / d,
194 (self.im*other.re - self.re*other.im) / d)
195
Alexander Belopolsky74135d02010-07-03 22:36:06 +0000196 def __rtruediv__(self, other):
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000197 other = ToComplex(other)
198 return other / self
199
200 def __pow__(self, n, z=None):
201 if z is not None:
Collin Winter6f2df4d2007-07-17 20:59:35 +0000202 raise TypeError('Complex does not support ternary pow()')
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000203 if IsComplex(n):
204 if n.im:
Collin Winter6f2df4d2007-07-17 20:59:35 +0000205 if self.im: raise TypeError('Complex to the Complex power')
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000206 else: return exp(math.log(self.re)*n)
207 n = n.re
208 r = pow(self.abs(), n)
209 phi = n*self.angle()
210 return Complex(math.cos(phi)*r, math.sin(phi)*r)
211
212 def __rpow__(self, base):
213 base = ToComplex(base)
214 return pow(base, self)
215
Guido van Rossum72ba6161996-07-30 19:02:01 +0000216def exp(z):
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000217 r = math.exp(z.re)
218 return Complex(math.cos(z.im)*r,math.sin(z.im)*r)
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000219
220
221def checkop(expr, a, b, value, fuzz = 1e-6):
Collin Winter6f2df4d2007-07-17 20:59:35 +0000222 print(' ', a, 'and', b, end=' ')
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000223 try:
224 result = eval(expr)
Alexander Belopolsky74135d02010-07-03 22:36:06 +0000225 except Exception as e:
226 print('!!\t!!\t!! error: {}'.format(e))
227 return
Collin Winter6f2df4d2007-07-17 20:59:35 +0000228 print('->', result)
Raymond Hettingerbdaad8c2005-04-09 14:55:07 +0000229 if isinstance(result, str) or isinstance(value, str):
230 ok = (result == value)
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000231 else:
232 ok = abs(result - value) <= fuzz
233 if not ok:
Collin Winter6f2df4d2007-07-17 20:59:35 +0000234 print('!!\t!!\t!! should be', value, 'diff', abs(result - value))
Guido van Rossume8769491992-08-13 12:14:11 +0000235
Guido van Rossume8769491992-08-13 12:14:11 +0000236def test():
Collin Winter6f2df4d2007-07-17 20:59:35 +0000237 print('test constructors')
Martin v. Löwis4a1e48c2005-04-09 10:51:19 +0000238 constructor_test = (
239 # "expect" is an array [re,im] "got" the Complex.
240 ( (0,0), Complex() ),
241 ( (0,0), Complex() ),
242 ( (1,0), Complex(1) ),
243 ( (0,1), Complex(0,1) ),
244 ( (1,2), Complex(Complex(1,2)) ),
245 ( (1,3), Complex(Complex(1,2),1) ),
246 ( (0,0), Complex(0,Complex(0,0)) ),
247 ( (3,4), Complex(3,Complex(4)) ),
248 ( (-1,3), Complex(1,Complex(3,2)) ),
249 ( (-7,6), Complex(Complex(1,2),Complex(4,8)) ) )
250 cnt = [0,0]
251 for t in constructor_test:
252 cnt[0] += 1
253 if ((t[0][0]!=t[1].re)or(t[0][1]!=t[1].im)):
Collin Winter6f2df4d2007-07-17 20:59:35 +0000254 print(" expected", t[0], "got", t[1])
Martin v. Löwis4a1e48c2005-04-09 10:51:19 +0000255 cnt[1] += 1
Collin Winter6f2df4d2007-07-17 20:59:35 +0000256 print(" ", cnt[1], "of", cnt[0], "tests failed")
Martin v. Löwis4a1e48c2005-04-09 10:51:19 +0000257 # test operators
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000258 testsuite = {
259 'a+b': [
260 (1, 10, 11),
261 (1, Complex(0,10), Complex(1,10)),
262 (Complex(0,10), 1, Complex(1,10)),
263 (Complex(0,10), Complex(1), Complex(1,10)),
264 (Complex(1), Complex(0,10), Complex(1,10)),
265 ],
266 'a-b': [
267 (1, 10, -9),
268 (1, Complex(0,10), Complex(1,-10)),
269 (Complex(0,10), 1, Complex(-1,10)),
270 (Complex(0,10), Complex(1), Complex(-1,10)),
271 (Complex(1), Complex(0,10), Complex(1,-10)),
272 ],
273 'a*b': [
274 (1, 10, 10),
275 (1, Complex(0,10), Complex(0, 10)),
276 (Complex(0,10), 1, Complex(0,10)),
277 (Complex(0,10), Complex(1), Complex(0,10)),
278 (Complex(1), Complex(0,10), Complex(0,10)),
279 ],
280 'a/b': [
281 (1., 10, 0.1),
282 (1, Complex(0,10), Complex(0, -0.1)),
283 (Complex(0, 10), 1, Complex(0, 10)),
284 (Complex(0, 10), Complex(1), Complex(0, 10)),
285 (Complex(1), Complex(0,10), Complex(0, -0.1)),
286 ],
287 'pow(a,b)': [
288 (1, 10, 1),
289 (1, Complex(0,10), 1),
290 (Complex(0,10), 1, Complex(0,10)),
291 (Complex(0,10), Complex(1), Complex(0,10)),
292 (Complex(1), Complex(0,10), 1),
293 (2, Complex(4,0), 16),
294 ],
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000295 }
Raymond Hettingerbdaad8c2005-04-09 14:55:07 +0000296 for expr in sorted(testsuite):
Collin Winter6f2df4d2007-07-17 20:59:35 +0000297 print(expr + ':')
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000298 t = (expr,)
299 for item in testsuite[expr]:
Raymond Hettingerbdaad8c2005-04-09 14:55:07 +0000300 checkop(*(t+item))
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000301
Guido van Rossume8769491992-08-13 12:14:11 +0000302
Guido van Rossum81a12bc1994-10-08 18:56:41 +0000303if __name__ == '__main__':
Andrew M. Kuchling946c53e2003-04-24 17:13:18 +0000304 test()