← Back to branch summary

~ubuntu-branches/ubuntu/maverick/python3.1/maverick

« back to all changes in this revision

Viewing changes to Objects/complexobject.c

  • Committer: Bazaar Package Importer
  • Author(s): Matthias Klose
  • Date: 2009-03-23 00:01:27 UTC
  • Revision ID: james.westby@ubuntu.com-20090323000127-5fstfxju4ufrhthq
Tags: upstream-3.1~a1+20090322
Import upstream version 3.1~a1+20090322

Show diffs side-by-side

added added

removed removed

Lines of Context:
Objects/complexobject.c
 
1
 
 
2
/* Complex object implementation */
 
3
 
 
4
/* Borrows heavily from floatobject.c */
 
5
 
 
6
/* Submitted by Jim Hugunin */
 
7
 
 
8
#include "Python.h"
 
9
#include "structmember.h"
 
10
 
 
11
#ifdef HAVE_IEEEFP_H
 
12
#include <ieeefp.h>
 
13
#endif
 
14
 
 
15
#ifndef WITHOUT_COMPLEX
 
16
 
 
17
/* Precisions used by repr() and str(), respectively.
 
18
 
 
19
   The repr() precision (17 significant decimal digits) is the minimal number
 
20
   that is guaranteed to have enough precision so that if the number is read
 
21
   back in the exact same binary value is recreated.  This is true for IEEE
 
22
   floating point by design, and also happens to work for all other modern
 
23
   hardware.
 
24
 
 
25
   The str() precision is chosen so that in most cases, the rounding noise
 
26
   created by various operations is suppressed, while giving plenty of
 
27
   precision for practical use.
 
28
*/
 
29
 
 
30
#define PREC_REPR       17
 
31
#define PREC_STR        12
 
32
 
 
33
/* elementary operations on complex numbers */
 
34
 
 
35
static Py_complex c_1 = {1., 0.};
 
36
 
 
37
Py_complex
 
38
c_sum(Py_complex a, Py_complex b)
 
39
{
 
40
        Py_complex r;
 
41
        r.real = a.real + b.real;
 
42
        r.imag = a.imag + b.imag;
 
43
        return r;
 
44
}
 
45
 
 
46
Py_complex
 
47
c_diff(Py_complex a, Py_complex b)
 
48
{
 
49
        Py_complex r;
 
50
        r.real = a.real - b.real;
 
51
        r.imag = a.imag - b.imag;
 
52
        return r;
 
53
}
 
54
 
 
55
Py_complex
 
56
c_neg(Py_complex a)
 
57
{
 
58
        Py_complex r;
 
59
        r.real = -a.real;
 
60
        r.imag = -a.imag;
 
61
        return r;
 
62
}
 
63
 
 
64
Py_complex
 
65
c_prod(Py_complex a, Py_complex b)
 
66
{
 
67
        Py_complex r;
 
68
        r.real = a.real*b.real - a.imag*b.imag;
 
69
        r.imag = a.real*b.imag + a.imag*b.real;
 
70
        return r;
 
71
}
 
72
 
 
73
Py_complex
 
74
c_quot(Py_complex a, Py_complex b)
 
75
{
 
76
        /******************************************************************
 
77
        This was the original algorithm.  It's grossly prone to spurious
 
78
        overflow and underflow errors.  It also merrily divides by 0 despite
 
79
        checking for that(!).  The code still serves a doc purpose here, as
 
80
        the algorithm following is a simple by-cases transformation of this
 
81
        one:
 
82
 
 
83
        Py_complex r;
 
84
        double d = b.real*b.real + b.imag*b.imag;
 
85
        if (d == 0.)
 
86
                errno = EDOM;
 
87
        r.real = (a.real*b.real + a.imag*b.imag)/d;
 
88
        r.imag = (a.imag*b.real - a.real*b.imag)/d;
 
89
        return r;
 
90
        ******************************************************************/
 
91
 
 
92
        /* This algorithm is better, and is pretty obvious:  first divide the
 
93
         * numerators and denominator by whichever of {b.real, b.imag} has
 
94
         * larger magnitude.  The earliest reference I found was to CACM
 
95
         * Algorithm 116 (Complex Division, Robert L. Smith, Stanford
 
96
         * University).  As usual, though, we're still ignoring all IEEE
 
97
         * endcases.
 
98
         */
 
99
         Py_complex r;  /* the result */
 
100
         const double abs_breal = b.real < 0 ? -b.real : b.real;
 
101
         const double abs_bimag = b.imag < 0 ? -b.imag : b.imag;
 
102
 
 
103
         if (abs_breal >= abs_bimag) {
 
104
                /* divide tops and bottom by b.real */
 
105
                if (abs_breal == 0.0) {
 
106
                        errno = EDOM;
 
107
                        r.real = r.imag = 0.0;
 
108
                }
 
109
                else {
 
110
                        const double ratio = b.imag / b.real;
 
111
                        const double denom = b.real + b.imag * ratio;
 
112
                        r.real = (a.real + a.imag * ratio) / denom;
 
113
                        r.imag = (a.imag - a.real * ratio) / denom;
 
114
                }
 
115
        }
 
116
        else {
 
117
                /* divide tops and bottom by b.imag */
 
118
                const double ratio = b.real / b.imag;
 
119
                const double denom = b.real * ratio + b.imag;
 
120
                assert(b.imag != 0.0);
 
121
                r.real = (a.real * ratio + a.imag) / denom;
 
122
                r.imag = (a.imag * ratio - a.real) / denom;
 
123
        }
 
124
        return r;
 
125
}
 
126
 
 
127
Py_complex
 
128
c_pow(Py_complex a, Py_complex b)
 
129
{
 
130
        Py_complex r;
 
131
        double vabs,len,at,phase;
 
132
        if (b.real == 0. && b.imag == 0.) {
 
133
                r.real = 1.;
 
134
                r.imag = 0.;
 
135
        }
 
136
        else if (a.real == 0. && a.imag == 0.) {
 
137
                if (b.imag != 0. || b.real < 0.)
 
138
                        errno = EDOM;
 
139
                r.real = 0.;
 
140
                r.imag = 0.;
 
141
        }
 
142
        else {
 
143
                vabs = hypot(a.real,a.imag);
 
144
                len = pow(vabs,b.real);
 
145
                at = atan2(a.imag, a.real);
 
146
                phase = at*b.real;
 
147
                if (b.imag != 0.0) {
 
148
                        len /= exp(at*b.imag);
 
149
                        phase += b.imag*log(vabs);
 
150
                }
 
151
                r.real = len*cos(phase);
 
152
                r.imag = len*sin(phase);
 
153
        }
 
154
        return r;
 
155
}
 
156
 
 
157
static Py_complex
 
158
c_powu(Py_complex x, long n)
 
159
{
 
160
        Py_complex r, p;
 
161
        long mask = 1;
 
162
        r = c_1;
 
163
        p = x;
 
164
        while (mask > 0 && n >= mask) {
 
165
                if (n & mask)
 
166
                        r = c_prod(r,p);
 
167
                mask <<= 1;
 
168
                p = c_prod(p,p);
 
169
        }
 
170
        return r;
 
171
}
 
172
 
 
173
static Py_complex
 
174
c_powi(Py_complex x, long n)
 
175
{
 
176
        Py_complex cn;
 
177
 
 
178
        if (n > 100 || n < -100) {
 
179
                cn.real = (double) n;
 
180
                cn.imag = 0.;
 
181
                return c_pow(x,cn);
 
182
        }
 
183
        else if (n > 0)
 
184
                return c_powu(x,n);
 
185
        else
 
186
                return c_quot(c_1,c_powu(x,-n));
 
187
 
 
188
}
 
189
 
 
190
double
 
191
c_abs(Py_complex z)
 
192
{
 
193
        /* sets errno = ERANGE on overflow;  otherwise errno = 0 */
 
194
        double result;
 
195
 
 
196
        if (!Py_IS_FINITE(z.real) || !Py_IS_FINITE(z.imag)) {
 
197
                /* C99 rules: if either the real or the imaginary part is an
 
198
                   infinity, return infinity, even if the other part is a
 
199
                   NaN. */
 
200
                if (Py_IS_INFINITY(z.real)) {
 
201
                        result = fabs(z.real);
 
202
                        errno = 0;
 
203
                        return result;
 
204
                }
 
205
                if (Py_IS_INFINITY(z.imag)) {
 
206
                        result = fabs(z.imag);
 
207
                        errno = 0;
 
208
                        return result;
 
209
                }
 
210
                /* either the real or imaginary part is a NaN,
 
211
                   and neither is infinite. Result should be NaN. */
 
212
                return Py_NAN;
 
213
        }
 
214
        result = hypot(z.real, z.imag);
 
215
        if (!Py_IS_FINITE(result))
 
216
                errno = ERANGE;
 
217
        else
 
218
                errno = 0;
 
219
        return result;
 
220
}
 
221
 
 
222
static PyObject *
 
223
complex_subtype_from_c_complex(PyTypeObject *type, Py_complex cval)
 
224
{
 
225
        PyObject *op;
 
226
 
 
227
        op = type->tp_alloc(type, 0);
 
228
        if (op != NULL)
 
229
                ((PyComplexObject *)op)->cval = cval;
 
230
        return op;
 
231
}
 
232
 
 
233
PyObject *
 
234
PyComplex_FromCComplex(Py_complex cval)
 
235
{
 
236
        register PyComplexObject *op;
 
237
 
 
238
        /* Inline PyObject_New */
 
239
        op = (PyComplexObject *) PyObject_MALLOC(sizeof(PyComplexObject));
 
240
        if (op == NULL)
 
241
                return PyErr_NoMemory();
 
242
        PyObject_INIT(op, &PyComplex_Type);
 
243
        op->cval = cval;
 
244
        return (PyObject *) op;
 
245
}
 
246
 
 
247
static PyObject *
 
248
complex_subtype_from_doubles(PyTypeObject *type, double real, double imag)
 
249
{
 
250
        Py_complex c;
 
251
        c.real = real;
 
252
        c.imag = imag;
 
253
        return complex_subtype_from_c_complex(type, c);
 
254
}
 
255
 
 
256
PyObject *
 
257
PyComplex_FromDoubles(double real, double imag)
 
258
{
 
259
        Py_complex c;
 
260
        c.real = real;
 
261
        c.imag = imag;
 
262
        return PyComplex_FromCComplex(c);
 
263
}
 
264
 
 
265
double
 
266
PyComplex_RealAsDouble(PyObject *op)
 
267
{
 
268
        if (PyComplex_Check(op)) {
 
269
                return ((PyComplexObject *)op)->cval.real;
 
270
        }
 
271
        else {
 
272
                return PyFloat_AsDouble(op);
 
273
        }
 
274
}
 
275
 
 
276
double
 
277
PyComplex_ImagAsDouble(PyObject *op)
 
278
{
 
279
        if (PyComplex_Check(op)) {
 
280
                return ((PyComplexObject *)op)->cval.imag;
 
281
        }
 
282
        else {
 
283
                return 0.0;
 
284
        }
 
285
}
 
286
 
 
287
Py_complex
 
288
PyComplex_AsCComplex(PyObject *op)
 
289
{
 
290
        Py_complex cv;
 
291
        PyObject *newop = NULL;
 
292
        static PyObject *complex_str = NULL;
 
293
 
 
294
        assert(op);
 
295
        /* If op is already of type PyComplex_Type, return its value */
 
296
        if (PyComplex_Check(op)) {
 
297
                return ((PyComplexObject *)op)->cval;
 
298
        }
 
299
        /* If not, use op's __complex__  method, if it exists */
 
300
        
 
301
        /* return -1 on failure */
 
302
        cv.real = -1.;
 
303
        cv.imag = 0.;
 
304
                
 
305
        if (complex_str == NULL) {
 
306
                if (!(complex_str = PyUnicode_FromString("__complex__")))
 
307
                        return cv;
 
308
        }
 
309
 
 
310
        {
 
311
                PyObject *complexfunc;
 
312
                complexfunc = _PyType_Lookup(op->ob_type, complex_str);
 
313
                /* complexfunc is a borrowed reference */
 
314
                if (complexfunc) {
 
315
                        newop = PyObject_CallFunctionObjArgs(complexfunc, op, NULL);
 
316
                        if (!newop)
 
317
                                return cv;
 
318
                }
 
319
        }
 
320
 
 
321
        if (newop) {
 
322
                if (!PyComplex_Check(newop)) {
 
323
                        PyErr_SetString(PyExc_TypeError,
 
324
                                "__complex__ should return a complex object");
 
325
                        Py_DECREF(newop);
 
326
                        return cv;
 
327
                }
 
328
                cv = ((PyComplexObject *)newop)->cval;
 
329
                Py_DECREF(newop);
 
330
                return cv;
 
331
        }
 
332
        /* If neither of the above works, interpret op as a float giving the
 
333
           real part of the result, and fill in the imaginary part as 0. */
 
334
        else {
 
335
                /* PyFloat_AsDouble will return -1 on failure */
 
336
                cv.real = PyFloat_AsDouble(op);
 
337
                return cv;
 
338
        }
 
339
}
 
340
 
 
341
static void
 
342
complex_dealloc(PyObject *op)
 
343
{
 
344
        op->ob_type->tp_free(op);
 
345
}
 
346
 
 
347
 
 
348
static void
 
349
complex_to_buf(char *buf, int bufsz, PyComplexObject *v, int precision)
 
350
{
 
351
        char format[32];
 
352
        if (v->cval.real == 0.) {
 
353
                if (!Py_IS_FINITE(v->cval.imag)) {
 
354
                        if (Py_IS_NAN(v->cval.imag))
 
355
                                strncpy(buf, "nan*j", 6);
 
356
                        else if (copysign(1, v->cval.imag) == 1)
 
357
                                strncpy(buf, "inf*j", 6);
 
358
                        else
 
359
                                strncpy(buf, "-inf*j", 7);
 
360
                }
 
361
                else {
 
362
                        PyOS_snprintf(format, sizeof(format), "%%.%ig", precision);
 
363
                        PyOS_ascii_formatd(buf, bufsz - 1, format, v->cval.imag);
 
364
                        strncat(buf, "j", 1);
 
365
                }
 
366
        } else {
 
367
                char re[64], im[64];
 
368
                /* Format imaginary part with sign, real part without */
 
369
                if (!Py_IS_FINITE(v->cval.real)) {
 
370
                        if (Py_IS_NAN(v->cval.real))
 
371
                                strncpy(re, "nan", 4);
 
372
                        /* else if (copysign(1, v->cval.real) == 1) */
 
373
                        else if (v->cval.real > 0)
 
374
                                strncpy(re, "inf", 4);
 
375
                        else
 
376
                                strncpy(re, "-inf", 5);
 
377
                }
 
378
                else {
 
379
                        PyOS_snprintf(format, sizeof(format), "%%.%ig", precision);
 
380
                        PyOS_ascii_formatd(re, sizeof(re), format, v->cval.real);
 
381
                }
 
382
                if (!Py_IS_FINITE(v->cval.imag)) {
 
383
                        if (Py_IS_NAN(v->cval.imag))
 
384
                                strncpy(im, "+nan*", 6);
 
385
                        /* else if (copysign(1, v->cval.imag) == 1) */
 
386
                        else if (v->cval.imag > 0)
 
387
                                strncpy(im, "+inf*", 6);
 
388
                        else
 
389
                                strncpy(im, "-inf*", 6);
 
390
                }
 
391
                else {
 
392
                        PyOS_snprintf(format, sizeof(format), "%%+.%ig", precision);
 
393
                        PyOS_ascii_formatd(im, sizeof(im), format, v->cval.imag);
 
394
                }
 
395
                PyOS_snprintf(buf, bufsz, "(%s%sj)", re, im);
 
396
        }
 
397
}
 
398
 
 
399
static PyObject *
 
400
complex_repr(PyComplexObject *v)
 
401
{
 
402
        char buf[100];
 
403
        complex_to_buf(buf, sizeof(buf), v, PREC_REPR);
 
404
        return PyUnicode_FromString(buf);
 
405
}
 
406
 
 
407
static PyObject *
 
408
complex_str(PyComplexObject *v)
 
409
{
 
410
        char buf[100];
 
411
        complex_to_buf(buf, sizeof(buf), v, PREC_STR);
 
412
        return PyUnicode_FromString(buf);
 
413
}
 
414
 
 
415
static long
 
416
complex_hash(PyComplexObject *v)
 
417
{
 
418
        long hashreal, hashimag, combined;
 
419
        hashreal = _Py_HashDouble(v->cval.real);
 
420
        if (hashreal == -1)
 
421
                return -1;
 
422
        hashimag = _Py_HashDouble(v->cval.imag);
 
423
        if (hashimag == -1)
 
424
                return -1;
 
425
        /* Note:  if the imaginary part is 0, hashimag is 0 now,
 
426
         * so the following returns hashreal unchanged.  This is
 
427
         * important because numbers of different types that
 
428
         * compare equal must have the same hash value, so that
 
429
         * hash(x + 0*j) must equal hash(x).
 
430
         */
 
431
        combined = hashreal + 1000003 * hashimag;
 
432
        if (combined == -1)
 
433
                combined = -2;
 
434
        return combined;
 
435
}
 
436
 
 
437
/* This macro may return! */
 
438
#define TO_COMPLEX(obj, c) \
 
439
        if (PyComplex_Check(obj)) \
 
440
                c = ((PyComplexObject *)(obj))->cval; \
 
441
        else if (to_complex(&(obj), &(c)) < 0) \
 
442
                return (obj)
 
443
 
 
444
static int
 
445
to_complex(PyObject **pobj, Py_complex *pc)
 
446
{
 
447
        PyObject *obj = *pobj;
 
448
 
 
449
        pc->real = pc->imag = 0.0;
 
450
        if (PyLong_Check(obj)) {
 
451
                pc->real = PyLong_AsDouble(obj);
 
452
                if (pc->real == -1.0 && PyErr_Occurred()) {
 
453
                        *pobj = NULL;
 
454
                        return -1;
 
455
                }
 
456
                return 0;
 
457
        }
 
458
        if (PyFloat_Check(obj)) {
 
459
                pc->real = PyFloat_AsDouble(obj);
 
460
                return 0;
 
461
        }
 
462
        Py_INCREF(Py_NotImplemented);
 
463
        *pobj = Py_NotImplemented;
 
464
        return -1;
 
465
}
 
466
                
 
467
 
 
468
static PyObject *
 
469
complex_add(PyObject *v, PyObject *w)
 
470
{
 
471
        Py_complex result;
 
472
        Py_complex a, b;
 
473
        TO_COMPLEX(v, a);
 
474
        TO_COMPLEX(w, b);
 
475
        PyFPE_START_PROTECT("complex_add", return 0)
 
476
        result = c_sum(a, b);
 
477
        PyFPE_END_PROTECT(result)
 
478
        return PyComplex_FromCComplex(result);
 
479
}
 
480
 
 
481
static PyObject *
 
482
complex_sub(PyObject *v, PyObject *w)
 
483
{
 
484
        Py_complex result;
 
485
        Py_complex a, b;
 
486
        TO_COMPLEX(v, a);
 
487
        TO_COMPLEX(w, b);
 
488
        PyFPE_START_PROTECT("complex_sub", return 0)
 
489
        result = c_diff(a, b);
 
490
        PyFPE_END_PROTECT(result)
 
491
        return PyComplex_FromCComplex(result);
 
492
}
 
493
 
 
494
static PyObject *
 
495
complex_mul(PyObject *v, PyObject *w)
 
496
{
 
497
        Py_complex result;
 
498
        Py_complex a, b;
 
499
        TO_COMPLEX(v, a);
 
500
        TO_COMPLEX(w, b);
 
501
        PyFPE_START_PROTECT("complex_mul", return 0)
 
502
        result = c_prod(a, b);
 
503
        PyFPE_END_PROTECT(result)
 
504
        return PyComplex_FromCComplex(result);
 
505
}
 
506
 
 
507
static PyObject *
 
508
complex_div(PyObject *v, PyObject *w)
 
509
{
 
510
        Py_complex quot;
 
511
        Py_complex a, b;
 
512
        TO_COMPLEX(v, a);
 
513
        TO_COMPLEX(w, b);
 
514
        PyFPE_START_PROTECT("complex_div", return 0)
 
515
        errno = 0;
 
516
        quot = c_quot(a, b);
 
517
        PyFPE_END_PROTECT(quot)
 
518
        if (errno == EDOM) {
 
519
                PyErr_SetString(PyExc_ZeroDivisionError, "complex division");
 
520
                return NULL;
 
521
        }
 
522
        return PyComplex_FromCComplex(quot);
 
523
}
 
524
 
 
525
static PyObject *
 
526
complex_remainder(PyObject *v, PyObject *w)
 
527
{
 
528
        PyErr_SetString(PyExc_TypeError,
 
529
                        "can't mod complex numbers.");
 
530
        return NULL;
 
531
}
 
532
 
 
533
 
 
534
static PyObject *
 
535
complex_divmod(PyObject *v, PyObject *w)
 
536
{
 
537
        PyErr_SetString(PyExc_TypeError,
 
538
                        "can't take floor or mod of complex number.");
 
539
        return NULL;
 
540
}
 
541
 
 
542
static PyObject *
 
543
complex_pow(PyObject *v, PyObject *w, PyObject *z)
 
544
{
 
545
        Py_complex p;
 
546
        Py_complex exponent;
 
547
        long int_exponent;
 
548
        Py_complex a, b;
 
549
        TO_COMPLEX(v, a);
 
550
        TO_COMPLEX(w, b);
 
551
 
 
552
        if (z != Py_None) {
 
553
                PyErr_SetString(PyExc_ValueError, "complex modulo");
 
554
                return NULL;
 
555
        }
 
556
        PyFPE_START_PROTECT("complex_pow", return 0)
 
557
        errno = 0;
 
558
        exponent = b;
 
559
        int_exponent = (long)exponent.real;
 
560
        if (exponent.imag == 0. && exponent.real == int_exponent)
 
561
                p = c_powi(a, int_exponent);
 
562
        else
 
563
                p = c_pow(a, exponent);
 
564
 
 
565
        PyFPE_END_PROTECT(p)
 
566
        Py_ADJUST_ERANGE2(p.real, p.imag);
 
567
        if (errno == EDOM) {
 
568
                PyErr_SetString(PyExc_ZeroDivisionError,
 
569
                                "0.0 to a negative or complex power");
 
570
                return NULL;
 
571
        }
 
572
        else if (errno == ERANGE) {
 
573
                PyErr_SetString(PyExc_OverflowError,
 
574
                                "complex exponentiation");
 
575
                return NULL;
 
576
        }
 
577
        return PyComplex_FromCComplex(p);
 
578
}
 
579
 
 
580
static PyObject *
 
581
complex_int_div(PyObject *v, PyObject *w)
 
582
{
 
583
        PyErr_SetString(PyExc_TypeError,
 
584
                        "can't take floor of complex number.");
 
585
        return NULL;
 
586
}
 
587
 
 
588
static PyObject *
 
589
complex_neg(PyComplexObject *v)
 
590
{
 
591
        Py_complex neg;
 
592
        neg.real = -v->cval.real;
 
593
        neg.imag = -v->cval.imag;
 
594
        return PyComplex_FromCComplex(neg);
 
595
}
 
596
 
 
597
static PyObject *
 
598
complex_pos(PyComplexObject *v)
 
599
{
 
600
        if (PyComplex_CheckExact(v)) {
 
601
                Py_INCREF(v);
 
602
                return (PyObject *)v;
 
603
        }
 
604
        else
 
605
                return PyComplex_FromCComplex(v->cval);
 
606
}
 
607
 
 
608
static PyObject *
 
609
complex_abs(PyComplexObject *v)
 
610
{
 
611
        double result;
 
612
 
 
613
        PyFPE_START_PROTECT("complex_abs", return 0)
 
614
        result = c_abs(v->cval);
 
615
        PyFPE_END_PROTECT(result)
 
616
 
 
617
        if (errno == ERANGE) {
 
618
                PyErr_SetString(PyExc_OverflowError,
 
619
                                "absolute value too large");
 
620
                return NULL;
 
621
        }
 
622
        return PyFloat_FromDouble(result);
 
623
}
 
624
 
 
625
static int
 
626
complex_bool(PyComplexObject *v)
 
627
{
 
628
        return v->cval.real != 0.0 || v->cval.imag != 0.0;
 
629
}
 
630
 
 
631
static PyObject *
 
632
complex_richcompare(PyObject *v, PyObject *w, int op)
 
633
{
 
634
        PyObject *res;
 
635
        Py_complex i, j;
 
636
        TO_COMPLEX(v, i);
 
637
        TO_COMPLEX(w, j);
 
638
 
 
639
        if (op != Py_EQ && op != Py_NE) {
 
640
                /* XXX Should eventually return NotImplemented */
 
641
                PyErr_SetString(PyExc_TypeError,
 
642
                        "no ordering relation is defined for complex numbers");
 
643
                return NULL;
 
644
        }
 
645
 
 
646
        if ((i.real == j.real && i.imag == j.imag) == (op == Py_EQ))
 
647
                res = Py_True;
 
648
        else
 
649
                res = Py_False;
 
650
 
 
651
        Py_INCREF(res);
 
652
        return res;
 
653
}
 
654
 
 
655
static PyObject *
 
656
complex_int(PyObject *v)
 
657
{
 
658
        PyErr_SetString(PyExc_TypeError,
 
659
                   "can't convert complex to int; use int(abs(z))");
 
660
        return NULL;
 
661
}
 
662
 
 
663
static PyObject *
 
664
complex_float(PyObject *v)
 
665
{
 
666
        PyErr_SetString(PyExc_TypeError,
 
667
                   "can't convert complex to float; use abs(z)");
 
668
        return NULL;
 
669
}
 
670
 
 
671
static PyObject *
 
672
complex_conjugate(PyObject *self)
 
673
{
 
674
        Py_complex c;
 
675
        c = ((PyComplexObject *)self)->cval;
 
676
        c.imag = -c.imag;
 
677
        return PyComplex_FromCComplex(c);
 
678
}
 
679
 
 
680
PyDoc_STRVAR(complex_conjugate_doc,
 
681
"complex.conjugate() -> complex\n"
 
682
"\n"
 
683
"Returns the complex conjugate of its argument. (3-4j).conjugate() == 3+4j.");
 
684
 
 
685
static PyObject *
 
686
complex_getnewargs(PyComplexObject *v)
 
687
{
 
688
        Py_complex c = v->cval;
 
689
        return Py_BuildValue("(dd)", c.real, c.imag);
 
690
}
 
691
 
 
692
#if 0
 
693
static PyObject *
 
694
complex_is_finite(PyObject *self)
 
695
{
 
696
        Py_complex c;
 
697
        c = ((PyComplexObject *)self)->cval;
 
698
        return PyBool_FromLong((long)(Py_IS_FINITE(c.real) &&
 
699
                                      Py_IS_FINITE(c.imag)));
 
700
}
 
701
 
 
702
PyDoc_STRVAR(complex_is_finite_doc,
 
703
"complex.is_finite() -> bool\n"
 
704
"\n"
 
705
"Returns True if the real and the imaginary part is finite.");
 
706
#endif
 
707
 
 
708
static PyMethodDef complex_methods[] = {
 
709
        {"conjugate",   (PyCFunction)complex_conjugate, METH_NOARGS,
 
710
         complex_conjugate_doc},
 
711
#if 0
 
712
        {"is_finite",   (PyCFunction)complex_is_finite, METH_NOARGS,
 
713
         complex_is_finite_doc},
 
714
#endif
 
715
        {"__getnewargs__",      (PyCFunction)complex_getnewargs,        METH_NOARGS},
 
716
        {NULL,          NULL}           /* sentinel */
 
717
};
 
718
 
 
719
static PyMemberDef complex_members[] = {
 
720
        {"real", T_DOUBLE, offsetof(PyComplexObject, cval.real), READONLY,
 
721
         "the real part of a complex number"},
 
722
        {"imag", T_DOUBLE, offsetof(PyComplexObject, cval.imag), READONLY,
 
723
         "the imaginary part of a complex number"},
 
724
        {0},
 
725
};
 
726
 
 
727
static PyObject *
 
728
complex_subtype_from_string(PyTypeObject *type, PyObject *v)
 
729
{
 
730
        const char *s, *start;
 
731
        char *end;
 
732
        double x=0.0, y=0.0, z;
 
733
        int got_re=0, got_im=0, got_bracket=0, done=0;
 
734
        int digit_or_dot;
 
735
        int sw_error=0;
 
736
        int sign;
 
737
        char buffer[256]; /* For errors */
 
738
        char s_buffer[256];
 
739
        Py_ssize_t len;
 
740
 
 
741
        if (PyUnicode_Check(v)) {
 
742
                if (PyUnicode_GET_SIZE(v) >= (Py_ssize_t)sizeof(s_buffer)) {
 
743
                        PyErr_SetString(PyExc_ValueError,
 
744
                                 "complex() literal too large to convert");
 
745
                        return NULL;
 
746
                }
 
747
                if (PyUnicode_EncodeDecimal(PyUnicode_AS_UNICODE(v),
 
748
                                            PyUnicode_GET_SIZE(v),
 
749
                                            s_buffer,
 
750
                                            NULL))
 
751
                        return NULL;
 
752
                s = s_buffer;
 
753
                len = strlen(s);
 
754
        }
 
755
        else if (PyObject_AsCharBuffer(v, &s, &len)) {
 
756
                PyErr_SetString(PyExc_TypeError,
 
757
                                "complex() arg is not a string");
 
758
                return NULL;
 
759
        }
 
760
 
 
761
        /* position on first nonblank */
 
762
        start = s;
 
763
        while (*s && isspace(Py_CHARMASK(*s)))
 
764
                s++;
 
765
        if (s[0] == '\0') {
 
766
                PyErr_SetString(PyExc_ValueError,
 
767
                                "complex() arg is an empty string");
 
768
                return NULL;
 
769
        }
 
770
        if (s[0] == '(') {
 
771
                /* Skip over possible bracket from repr(). */
 
772
                got_bracket = 1;
 
773
                s++;
 
774
                while (*s && isspace(Py_CHARMASK(*s)))
 
775
                        s++;
 
776
        }
 
777
 
 
778
        z = -1.0;
 
779
        sign = 1;
 
780
        do {
 
781
 
 
782
                switch (*s) {
 
783
 
 
784
                case '\0':
 
785
                        if (s-start != len) {
 
786
                                PyErr_SetString(
 
787
                                        PyExc_ValueError,
 
788
                                        "complex() arg contains a null byte");
 
789
                                return NULL;
 
790
                        }
 
791
                        if(!done) sw_error=1;
 
792
                        break;
 
793
 
 
794
                case ')':
 
795
                        if (!got_bracket || !(got_re || got_im)) {
 
796
                                sw_error=1;
 
797
                                break;
 
798
                        }
 
799
                        got_bracket=0;
 
800
                        done=1;
 
801
                        s++;
 
802
                        while (*s && isspace(Py_CHARMASK(*s)))
 
803
                                s++;
 
804
                        if (*s) sw_error=1;
 
805
                        break;
 
806
 
 
807
                case '-':
 
808
                        sign = -1;
 
809
                                /* Fallthrough */
 
810
                case '+':
 
811
                        if (done)  sw_error=1;
 
812
                        s++;
 
813
                        if  (  *s=='\0'||*s=='+'||*s=='-'||*s==')'||
 
814
                               isspace(Py_CHARMASK(*s))  )  sw_error=1;
 
815
                        break;
 
816
 
 
817
                case 'J':
 
818
                case 'j':
 
819
                        if (got_im || done) {
 
820
                                sw_error = 1;
 
821
                                break;
 
822
                        }
 
823
                        if  (z<0.0) {
 
824
                                y=sign;
 
825
                        }
 
826
                        else{
 
827
                                y=sign*z;
 
828
                        }
 
829
                        got_im=1;
 
830
                        s++;
 
831
                        if  (*s!='+' && *s!='-' )
 
832
                                done=1;
 
833
                        break;
 
834
 
 
835
                default:
 
836
                        if (isspace(Py_CHARMASK(*s))) {
 
837
                                while (*s && isspace(Py_CHARMASK(*s)))
 
838
                                        s++;
 
839
                                if (*s && *s != ')')
 
840
                                        sw_error=1;
 
841
                                else
 
842
                                        done = 1;
 
843
                                break;
 
844
                        }
 
845
                        digit_or_dot =
 
846
                                (*s=='.' || isdigit(Py_CHARMASK(*s)));
 
847
                        if  (done||!digit_or_dot) {
 
848
                                sw_error=1;
 
849
                                break;
 
850
                        }
 
851
                        errno = 0;
 
852
                        PyFPE_START_PROTECT("strtod", return 0)
 
853
                                z = PyOS_ascii_strtod(s, &end) ;
 
854
                        PyFPE_END_PROTECT(z)
 
855
                                if (errno != 0) {
 
856
                                        PyOS_snprintf(buffer, sizeof(buffer),
 
857
                                          "float() out of range: %.150s", s);
 
858
                                        PyErr_SetString(
 
859
                                                PyExc_ValueError,
 
860
                                                buffer);
 
861
                                        return NULL;
 
862
                                }
 
863
                        s=end;
 
864
                        if  (*s=='J' || *s=='j') {
 
865
 
 
866
                                break;
 
867
                        }
 
868
                        if  (got_re) {
 
869
                                sw_error=1;
 
870
                                break;
 
871
                        }
 
872
 
 
873
                                /* accept a real part */
 
874
                        x=sign*z;
 
875
                        got_re=1;
 
876
                        if  (got_im)  done=1;
 
877
                        z = -1.0;
 
878
                        sign = 1;
 
879
                        break;
 
880
 
 
881
                }  /* end of switch  */
 
882
 
 
883
        } while (s - start < len && !sw_error);
 
884
 
 
885
        if (sw_error || got_bracket) {
 
886
                PyErr_SetString(PyExc_ValueError,
 
887
                                "complex() arg is a malformed string");
 
888
                return NULL;
 
889
        }
 
890
 
 
891
        return complex_subtype_from_doubles(type, x, y);
 
892
}
 
893
 
 
894
static PyObject *
 
895
complex_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
 
896
{
 
897
        PyObject *r, *i, *tmp, *f;
 
898
        PyNumberMethods *nbr, *nbi = NULL;
 
899
        Py_complex cr, ci;
 
900
        int own_r = 0;
 
901
        int cr_is_complex = 0;
 
902
        int ci_is_complex = 0;
 
903
        static PyObject *complexstr;
 
904
        static char *kwlist[] = {"real", "imag", 0};
 
905
 
 
906
        r = Py_False;
 
907
        i = NULL;
 
908
        if (!PyArg_ParseTupleAndKeywords(args, kwds, "|OO:complex", kwlist,
 
909
                                         &r, &i))
 
910
                return NULL;
 
911
 
 
912
        /* Special-case for a single argument when type(arg) is complex. */
 
913
        if (PyComplex_CheckExact(r) && i == NULL &&
 
914
            type == &PyComplex_Type) {
 
915
                /* Note that we can't know whether it's safe to return
 
916
                   a complex *subclass* instance as-is, hence the restriction
 
917
                   to exact complexes here.  If either the input or the
 
918
                   output is a complex subclass, it will be handled below 
 
919
                   as a non-orthogonal vector.  */
 
920
                Py_INCREF(r);
 
921
                return r;
 
922
        }
 
923
        if (PyUnicode_Check(r)) {
 
924
                if (i != NULL) {
 
925
                        PyErr_SetString(PyExc_TypeError,
 
926
                                        "complex() can't take second arg"
 
927
                                        " if first is a string");
 
928
                        return NULL;
 
929
                }
 
930
                return complex_subtype_from_string(type, r);
 
931
        }
 
932
        if (i != NULL && PyUnicode_Check(i)) {
 
933
                PyErr_SetString(PyExc_TypeError,
 
934
                                "complex() second arg can't be a string");
 
935
                return NULL;
 
936
        }
 
937
 
 
938
        /* XXX Hack to support classes with __complex__ method */
 
939
        if (complexstr == NULL) {
 
940
                complexstr = PyUnicode_InternFromString("__complex__");
 
941
                if (complexstr == NULL)
 
942
                        return NULL;
 
943
        }
 
944
        f = PyObject_GetAttr(r, complexstr);
 
945
        if (f == NULL)
 
946
                PyErr_Clear();
 
947
        else {
 
948
                PyObject *args = PyTuple_New(0);
 
949
                if (args == NULL)
 
950
                        return NULL;
 
951
                r = PyEval_CallObject(f, args);
 
952
                Py_DECREF(args);
 
953
                Py_DECREF(f);
 
954
                if (r == NULL)
 
955
                        return NULL;
 
956
                own_r = 1;
 
957
        }
 
958
        nbr = r->ob_type->tp_as_number;
 
959
        if (i != NULL)
 
960
                nbi = i->ob_type->tp_as_number;
 
961
        if (nbr == NULL || nbr->nb_float == NULL ||
 
962
            ((i != NULL) && (nbi == NULL || nbi->nb_float == NULL))) {
 
963
                PyErr_SetString(PyExc_TypeError,
 
964
                           "complex() argument must be a string or a number");
 
965
                if (own_r) {
 
966
                        Py_DECREF(r);
 
967
                }
 
968
                return NULL;
 
969
        }
 
970
 
 
971
        /* If we get this far, then the "real" and "imag" parts should
 
972
           both be treated as numbers, and the constructor should return a
 
973
           complex number equal to (real + imag*1j).
 
974
 
 
975
           Note that we do NOT assume the input to already be in canonical
 
976
           form; the "real" and "imag" parts might themselves be complex
 
977
           numbers, which slightly complicates the code below. */
 
978
        if (PyComplex_Check(r)) {
 
979
                /* Note that if r is of a complex subtype, we're only
 
980
                   retaining its real & imag parts here, and the return
 
981
                   value is (properly) of the builtin complex type. */
 
982
                cr = ((PyComplexObject*)r)->cval;
 
983
                cr_is_complex = 1;
 
984
                if (own_r) {
 
985
                        Py_DECREF(r);
 
986
                }
 
987
        }
 
988
        else {
 
989
                /* The "real" part really is entirely real, and contributes
 
990
                   nothing in the imaginary direction.  
 
991
                   Just treat it as a double. */
 
992
                tmp = PyNumber_Float(r);
 
993
                if (own_r) {
 
994
                        /* r was a newly created complex number, rather
 
995
                           than the original "real" argument. */
 
996
                        Py_DECREF(r);
 
997
                }
 
998
                if (tmp == NULL)
 
999
                        return NULL;
 
1000
                if (!PyFloat_Check(tmp)) {
 
1001
                        PyErr_SetString(PyExc_TypeError,
 
1002
                                        "float(r) didn't return a float");
 
1003
                        Py_DECREF(tmp);
 
1004
                        return NULL;
 
1005
                }
 
1006
                cr.real = PyFloat_AsDouble(tmp);
 
1007
                cr.imag = 0.0; /* Shut up compiler warning */
 
1008
                Py_DECREF(tmp);
 
1009
        }
 
1010
        if (i == NULL) {
 
1011
                ci.real = 0.0;
 
1012
        }
 
1013
        else if (PyComplex_Check(i)) {
 
1014
                ci = ((PyComplexObject*)i)->cval;
 
1015
                ci_is_complex = 1;
 
1016
        } else {
 
1017
                /* The "imag" part really is entirely imaginary, and
 
1018
                   contributes nothing in the real direction.
 
1019
                   Just treat it as a double. */
 
1020
                tmp = (*nbi->nb_float)(i);
 
1021
                if (tmp == NULL)
 
1022
                        return NULL;
 
1023
                ci.real = PyFloat_AsDouble(tmp);
 
1024
                Py_DECREF(tmp);
 
1025
        }
 
1026
        /*  If the input was in canonical form, then the "real" and "imag"
 
1027
            parts are real numbers, so that ci.imag and cr.imag are zero.
 
1028
            We need this correction in case they were not real numbers. */
 
1029
 
 
1030
        if (ci_is_complex) {
 
1031
                cr.real -= ci.imag;
 
1032
        }
 
1033
        if (cr_is_complex) {
 
1034
                ci.real += cr.imag;
 
1035
        }
 
1036
        return complex_subtype_from_doubles(type, cr.real, ci.real);
 
1037
}
 
1038
 
 
1039
PyDoc_STRVAR(complex_doc,
 
1040
"complex(real[, imag]) -> complex number\n"
 
1041
"\n"
 
1042
"Create a complex number from a real part and an optional imaginary part.\n"
 
1043
"This is equivalent to (real + imag*1j) where imag defaults to 0.");
 
1044
 
 
1045
static PyNumberMethods complex_as_number = {
 
1046
        (binaryfunc)complex_add,                /* nb_add */
 
1047
        (binaryfunc)complex_sub,                /* nb_subtract */
 
1048
        (binaryfunc)complex_mul,                /* nb_multiply */
 
1049
        (binaryfunc)complex_remainder,          /* nb_remainder */
 
1050
        (binaryfunc)complex_divmod,             /* nb_divmod */
 
1051
        (ternaryfunc)complex_pow,               /* nb_power */
 
1052
        (unaryfunc)complex_neg,                 /* nb_negative */
 
1053
        (unaryfunc)complex_pos,                 /* nb_positive */
 
1054
        (unaryfunc)complex_abs,                 /* nb_absolute */
 
1055
        (inquiry)complex_bool,                  /* nb_bool */
 
1056
        0,                                      /* nb_invert */
 
1057
        0,                                      /* nb_lshift */
 
1058
        0,                                      /* nb_rshift */
 
1059
        0,                                      /* nb_and */
 
1060
        0,                                      /* nb_xor */
 
1061
        0,                                      /* nb_or */
 
1062
        complex_int,                            /* nb_int */
 
1063
        0,                                      /* nb_reserved */
 
1064
        complex_float,                          /* nb_float */
 
1065
        0,                                      /* nb_inplace_add */
 
1066
        0,                                      /* nb_inplace_subtract */
 
1067
        0,                                      /* nb_inplace_multiply*/
 
1068
        0,                                      /* nb_inplace_remainder */
 
1069
        0,                                      /* nb_inplace_power */
 
1070
        0,                                      /* nb_inplace_lshift */
 
1071
        0,                                      /* nb_inplace_rshift */
 
1072
        0,                                      /* nb_inplace_and */
 
1073
        0,                                      /* nb_inplace_xor */
 
1074
        0,                                      /* nb_inplace_or */
 
1075
        (binaryfunc)complex_int_div,            /* nb_floor_divide */
 
1076
        (binaryfunc)complex_div,                /* nb_true_divide */
 
1077
        0,                                      /* nb_inplace_floor_divide */
 
1078
        0,                                      /* nb_inplace_true_divide */
 
1079
};
 
1080
 
 
1081
PyTypeObject PyComplex_Type = {
 
1082
        PyVarObject_HEAD_INIT(&PyType_Type, 0)
 
1083
        "complex",
 
1084
        sizeof(PyComplexObject),
 
1085
        0,
 
1086
        complex_dealloc,                        /* tp_dealloc */
 
1087
        0,                                      /* tp_print */
 
1088
        0,                                      /* tp_getattr */
 
1089
        0,                                      /* tp_setattr */
 
1090
        0,                                      /* tp_reserved */
 
1091
        (reprfunc)complex_repr,                 /* tp_repr */
 
1092
        &complex_as_number,                     /* tp_as_number */
 
1093
        0,                                      /* tp_as_sequence */
 
1094
        0,                                      /* tp_as_mapping */
 
1095
        (hashfunc)complex_hash,                 /* tp_hash */
 
1096
        0,                                      /* tp_call */
 
1097
        (reprfunc)complex_str,                  /* tp_str */
 
1098
        PyObject_GenericGetAttr,                /* tp_getattro */
 
1099
        0,                                      /* tp_setattro */
 
1100
        0,                                      /* tp_as_buffer */
 
1101
        Py_TPFLAGS_DEFAULT | Py_TPFLAGS_BASETYPE, /* tp_flags */
 
1102
        complex_doc,                            /* tp_doc */
 
1103
        0,                                      /* tp_traverse */
 
1104
        0,                                      /* tp_clear */
 
1105
        complex_richcompare,                    /* tp_richcompare */
 
1106
        0,                                      /* tp_weaklistoffset */
 
1107
        0,                                      /* tp_iter */
 
1108
        0,                                      /* tp_iternext */
 
1109
        complex_methods,                        /* tp_methods */
 
1110
        complex_members,                        /* tp_members */
 
1111
        0,                                      /* tp_getset */
 
1112
        0,                                      /* tp_base */
 
1113
        0,                                      /* tp_dict */
 
1114
        0,                                      /* tp_descr_get */
 
1115
        0,                                      /* tp_descr_set */
 
1116
        0,                                      /* tp_dictoffset */
 
1117
        0,                                      /* tp_init */
 
1118
        PyType_GenericAlloc,                    /* tp_alloc */
 
1119
        complex_new,                            /* tp_new */
 
1120
        PyObject_Del,                           /* tp_free */
 
1121
};
 
1122
 
 
1123
#endif