LCOV - code coverage report
Current view: top level - py - objcomplex.c (source / functions) Hit Total Coverage
Test: unix_coverage_v1.23.0-410-g82e69df33.info Lines: 125 125 100.0 %
Date: 2024-10-10 13:49:24 Functions: 8 8 100.0 %
Branches: 71 75 94.7 %

           Branch data     Line data    Source code
       1                 :            : /*
       2                 :            :  * This file is part of the MicroPython project, http://micropython.org/
       3                 :            :  *
       4                 :            :  * The MIT License (MIT)
       5                 :            :  *
       6                 :            :  * Copyright (c) 2013, 2014 Damien P. George
       7                 :            :  *
       8                 :            :  * Permission is hereby granted, free of charge, to any person obtaining a copy
       9                 :            :  * of this software and associated documentation files (the "Software"), to deal
      10                 :            :  * in the Software without restriction, including without limitation the rights
      11                 :            :  * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
      12                 :            :  * copies of the Software, and to permit persons to whom the Software is
      13                 :            :  * furnished to do so, subject to the following conditions:
      14                 :            :  *
      15                 :            :  * The above copyright notice and this permission notice shall be included in
      16                 :            :  * all copies or substantial portions of the Software.
      17                 :            :  *
      18                 :            :  * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
      19                 :            :  * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
      20                 :            :  * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
      21                 :            :  * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
      22                 :            :  * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
      23                 :            :  * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
      24                 :            :  * THE SOFTWARE.
      25                 :            :  */
      26                 :            : 
      27                 :            : #include <stdlib.h>
      28                 :            : #include <stdio.h>
      29                 :            : #include <assert.h>
      30                 :            : 
      31                 :            : #include "py/parsenum.h"
      32                 :            : #include "py/runtime.h"
      33                 :            : 
      34                 :            : #if MICROPY_PY_BUILTINS_COMPLEX
      35                 :            : 
      36                 :            : #include <math.h>
      37                 :            : #include "py/formatfloat.h"
      38                 :            : 
      39                 :            : typedef struct _mp_obj_complex_t {
      40                 :            :     mp_obj_base_t base;
      41                 :            :     mp_float_t real;
      42                 :            :     mp_float_t imag;
      43                 :            : } mp_obj_complex_t;
      44                 :            : 
      45                 :        344 : static void complex_print(const mp_print_t *print, mp_obj_t o_in, mp_print_kind_t kind) {
      46                 :        344 :     (void)kind;
      47                 :        344 :     mp_obj_complex_t *o = MP_OBJ_TO_PTR(o_in);
      48                 :            :     #if MICROPY_FLOAT_IMPL == MICROPY_FLOAT_IMPL_FLOAT
      49                 :            :     char buf[16];
      50                 :            :     #if MICROPY_OBJ_REPR == MICROPY_OBJ_REPR_C
      51                 :            :     const int precision = 6;
      52                 :            :     #else
      53                 :            :     const int precision = 7;
      54                 :            :     #endif
      55                 :            :     #else
      56                 :        344 :     char buf[32];
      57                 :        344 :     const int precision = 16;
      58                 :            :     #endif
      59         [ +  + ]:        344 :     if (o->real == 0) {
      60                 :         72 :         mp_format_float(o->imag, buf, sizeof(buf), 'g', precision, '\0');
      61                 :         72 :         mp_printf(print, "%sj", buf);
      62                 :            :     } else {
      63                 :        272 :         mp_format_float(o->real, buf, sizeof(buf), 'g', precision, '\0');
      64                 :        272 :         mp_printf(print, "(%s", buf);
      65   [ +  +  +  + ]:        272 :         if (o->imag >= 0 || isnan(o->imag)) {
      66                 :        172 :             mp_print_str(print, "+");
      67                 :            :         }
      68                 :        272 :         mp_format_float(o->imag, buf, sizeof(buf), 'g', precision, '\0');
      69                 :        272 :         mp_printf(print, "%sj)", buf);
      70                 :            :     }
      71                 :        344 : }
      72                 :            : 
      73                 :        508 : static mp_obj_t complex_make_new(const mp_obj_type_t *type_in, size_t n_args, size_t n_kw, const mp_obj_t *args) {
      74                 :        508 :     (void)type_in;
      75                 :        508 :     mp_arg_check_num(n_args, n_kw, 0, 2, false);
      76                 :            : 
      77      [ +  +  + ]:        508 :     switch (n_args) {
      78                 :          4 :         case 0:
      79                 :          4 :             return mp_obj_new_complex(0, 0);
      80                 :            : 
      81                 :        104 :         case 1:
      82   [ +  +  +  +  :        104 :             if (mp_obj_is_str(args[0])) {
                   +  + ]
      83                 :            :                 // a string, parse it
      84                 :         72 :                 size_t l;
      85                 :         72 :                 const char *s = mp_obj_str_get_data(args[0], &l);
      86                 :         72 :                 return mp_parse_num_complex(s, l, NULL);
      87   [ +  +  +  + ]:         32 :             } else if (mp_obj_is_type(args[0], &mp_type_complex)) {
      88                 :            :                 // a complex, just return it
      89                 :            :                 return args[0];
      90                 :            :             } else {
      91                 :         28 :                 mp_float_t real, imag;
      92                 :         28 :                 mp_obj_get_complex(args[0], &real, &imag);
      93                 :         16 :                 return mp_obj_new_complex(real, imag);
      94                 :            :             }
      95                 :            : 
      96                 :            :         case 2:
      97                 :            :         default: {
      98                 :        400 :             mp_float_t real, imag;
      99   [ +  +  +  + ]:        400 :             if (mp_obj_is_type(args[0], &mp_type_complex)) {
     100                 :          4 :                 mp_obj_complex_get(args[0], &real, &imag);
     101                 :            :             } else {
     102                 :        396 :                 real = mp_obj_get_float(args[0]);
     103                 :        396 :                 imag = 0;
     104                 :            :             }
     105   [ +  +  +  + ]:        400 :             if (mp_obj_is_type(args[1], &mp_type_complex)) {
     106                 :          4 :                 mp_float_t real2, imag2;
     107                 :          4 :                 mp_obj_complex_get(args[1], &real2, &imag2);
     108                 :          4 :                 real -= imag2;
     109                 :          4 :                 imag += real2;
     110                 :            :             } else {
     111                 :        396 :                 imag += mp_obj_get_float(args[1]);
     112                 :            :             }
     113                 :        400 :             return mp_obj_new_complex(real, imag);
     114                 :            :         }
     115                 :            :     }
     116                 :            : }
     117                 :            : 
     118                 :       2884 : static mp_obj_t complex_unary_op(mp_unary_op_t op, mp_obj_t o_in) {
     119                 :       2884 :     mp_obj_complex_t *o = MP_OBJ_TO_PTR(o_in);
     120   [ +  +  +  +  :       2884 :     switch (op) {
                   +  + ]
     121                 :          8 :         case MP_UNARY_OP_BOOL:
     122   [ +  -  +  + ]:          8 :             return mp_obj_new_bool(o->real != 0 || o->imag != 0);
     123                 :          8 :         case MP_UNARY_OP_HASH:
     124                 :          8 :             return MP_OBJ_NEW_SMALL_INT(mp_float_hash(o->real) ^ mp_float_hash(o->imag));
     125                 :            :         case MP_UNARY_OP_POSITIVE:
     126                 :            :             return o_in;
     127                 :          4 :         case MP_UNARY_OP_NEGATIVE:
     128                 :          4 :             return mp_obj_new_complex(-o->real, -o->imag);
     129                 :          8 :         case MP_UNARY_OP_ABS:
     130                 :          8 :             return mp_obj_new_float(MICROPY_FLOAT_C_FUN(sqrt)(o->real * o->real + o->imag * o->imag));
     131                 :       2852 :         default:
     132                 :       2852 :             return MP_OBJ_NULL;      // op not supported
     133                 :            :     }
     134                 :            : }
     135                 :            : 
     136                 :       1304 : static mp_obj_t complex_binary_op(mp_binary_op_t op, mp_obj_t lhs_in, mp_obj_t rhs_in) {
     137                 :       1304 :     mp_obj_complex_t *lhs = MP_OBJ_TO_PTR(lhs_in);
     138                 :       1304 :     return mp_obj_complex_binary_op(op, lhs->real, lhs->imag, rhs_in);
     139                 :            : }
     140                 :            : 
     141                 :       2484 : static void complex_attr(mp_obj_t self_in, qstr attr, mp_obj_t *dest) {
     142         [ +  + ]:       2484 :     if (dest[0] != MP_OBJ_NULL) {
     143                 :            :         // not load attribute
     144                 :            :         return;
     145                 :            :     }
     146                 :       2480 :     mp_obj_complex_t *self = MP_OBJ_TO_PTR(self_in);
     147         [ +  + ]:       2480 :     if (attr == MP_QSTR_real) {
     148                 :       1240 :         dest[0] = mp_obj_new_float(self->real);
     149         [ +  - ]:       1240 :     } else if (attr == MP_QSTR_imag) {
     150                 :       1240 :         dest[0] = mp_obj_new_float(self->imag);
     151                 :            :     }
     152                 :            : }
     153                 :            : 
     154                 :            : MP_DEFINE_CONST_OBJ_TYPE(
     155                 :            :     mp_type_complex, MP_QSTR_complex, MP_TYPE_FLAG_EQ_NOT_REFLEXIVE | MP_TYPE_FLAG_EQ_CHECKS_OTHER_TYPE,
     156                 :            :     make_new, complex_make_new,
     157                 :            :     print, complex_print,
     158                 :            :     unary_op, complex_unary_op,
     159                 :            :     binary_op, complex_binary_op,
     160                 :            :     attr, complex_attr
     161                 :            :     );
     162                 :            : 
     163                 :       1988 : mp_obj_t mp_obj_new_complex(mp_float_t real, mp_float_t imag) {
     164                 :       1988 :     mp_obj_complex_t *o = mp_obj_malloc(mp_obj_complex_t, &mp_type_complex);
     165                 :       1988 :     o->real = real;
     166                 :       1988 :     o->imag = imag;
     167                 :       1988 :     return MP_OBJ_FROM_PTR(o);
     168                 :            : }
     169                 :            : 
     170                 :       2864 : void mp_obj_complex_get(mp_obj_t self_in, mp_float_t *real, mp_float_t *imag) {
     171   [ +  -  -  + ]:       2864 :     assert(mp_obj_is_type(self_in, &mp_type_complex));
     172                 :       2864 :     mp_obj_complex_t *self = MP_OBJ_TO_PTR(self_in);
     173                 :       2864 :     *real = self->real;
     174                 :       2864 :     *imag = self->imag;
     175                 :       2864 : }
     176                 :            : 
     177                 :       1384 : mp_obj_t mp_obj_complex_binary_op(mp_binary_op_t op, mp_float_t lhs_real, mp_float_t lhs_imag, mp_obj_t rhs_in) {
     178                 :       1384 :     mp_float_t rhs_real, rhs_imag;
     179         [ +  + ]:       1384 :     if (!mp_obj_get_complex_maybe(rhs_in, &rhs_real, &rhs_imag)) {
     180                 :            :         return MP_OBJ_NULL; // op not supported
     181                 :            :     }
     182                 :            : 
     183   [ +  +  +  +  :       1376 :     switch (op) {
             +  +  +  + ]
     184                 :         72 :         case MP_BINARY_OP_ADD:
     185                 :            :         case MP_BINARY_OP_INPLACE_ADD:
     186                 :         72 :             lhs_real += rhs_real;
     187                 :         72 :             lhs_imag += rhs_imag;
     188                 :         72 :             break;
     189                 :          8 :         case MP_BINARY_OP_SUBTRACT:
     190                 :            :         case MP_BINARY_OP_INPLACE_SUBTRACT:
     191                 :          8 :             lhs_real -= rhs_real;
     192                 :          8 :             lhs_imag -= rhs_imag;
     193                 :          8 :             break;
     194                 :            :         case MP_BINARY_OP_MULTIPLY:
     195                 :            :         case MP_BINARY_OP_INPLACE_MULTIPLY: {
     196                 :         36 :             mp_float_t real;
     197                 :         36 :         multiply:
     198                 :         36 :             real = lhs_real * rhs_real - lhs_imag * rhs_imag;
     199                 :         36 :             lhs_imag = lhs_real * rhs_imag + lhs_imag * rhs_real;
     200                 :         36 :             lhs_real = real;
     201                 :         36 :             break;
     202                 :            :         }
     203                 :            :         case MP_BINARY_OP_FLOOR_DIVIDE:
     204                 :            :         case MP_BINARY_OP_INPLACE_FLOOR_DIVIDE:
     205                 :          4 :             mp_raise_TypeError(MP_ERROR_TEXT("can't truncate-divide a complex number"));
     206                 :            : 
     207                 :         16 :         case MP_BINARY_OP_TRUE_DIVIDE:
     208                 :            :         case MP_BINARY_OP_INPLACE_TRUE_DIVIDE:
     209         [ +  + ]:         16 :             if (rhs_imag == 0) {
     210         [ +  + ]:          8 :                 if (rhs_real == 0) {
     211                 :          4 :                     mp_raise_msg(&mp_type_ZeroDivisionError, MP_ERROR_TEXT("complex divide by zero"));
     212                 :            :                 }
     213                 :          4 :                 lhs_real /= rhs_real;
     214                 :          4 :                 lhs_imag /= rhs_real;
     215         [ +  + ]:          8 :             } else if (rhs_real == 0) {
     216                 :          4 :                 mp_float_t real = lhs_imag / rhs_imag;
     217                 :          4 :                 lhs_imag = -lhs_real / rhs_imag;
     218                 :          4 :                 lhs_real = real;
     219                 :            :             } else {
     220                 :          4 :                 mp_float_t rhs_len_sq = rhs_real * rhs_real + rhs_imag * rhs_imag;
     221                 :          4 :                 rhs_real /= rhs_len_sq;
     222                 :          4 :                 rhs_imag /= -rhs_len_sq;
     223                 :          4 :                 goto multiply;
     224                 :            :             }
     225                 :            :             break;
     226                 :            : 
     227                 :         36 :         case MP_BINARY_OP_POWER:
     228                 :            :         case MP_BINARY_OP_INPLACE_POWER: {
     229                 :            :             // z1**z2 = exp(z2*ln(z1))
     230                 :            :             //        = exp(z2*(ln(|z1|)+i*arg(z1)))
     231                 :            :             //        = exp( (x2*ln1 - y2*arg1) + i*(y2*ln1 + x2*arg1) )
     232                 :            :             //        = exp(x3 + i*y3)
     233                 :            :             //        = exp(x3)*(cos(y3) + i*sin(y3))
     234                 :         36 :             mp_float_t abs1 = MICROPY_FLOAT_C_FUN(sqrt)(lhs_real * lhs_real + lhs_imag * lhs_imag);
     235         [ +  + ]:         36 :             if (abs1 == 0) {
     236   [ +  +  +  + ]:         20 :                 if (rhs_imag == 0 && rhs_real >= 0) {
     237         [ +  + ]:         12 :                     lhs_real = (rhs_real == 0);
     238                 :            :                 } else {
     239                 :          8 :                     mp_raise_msg(&mp_type_ZeroDivisionError, MP_ERROR_TEXT("0.0 to a complex power"));
     240                 :            :                 }
     241                 :            :             } else {
     242                 :         16 :                 mp_float_t ln1 = MICROPY_FLOAT_C_FUN(log)(abs1);
     243                 :         16 :                 mp_float_t arg1 = MICROPY_FLOAT_C_FUN(atan2)(lhs_imag, lhs_real);
     244                 :         16 :                 mp_float_t x3 = rhs_real * ln1 - rhs_imag * arg1;
     245                 :         16 :                 mp_float_t y3 = rhs_imag * ln1 + rhs_real * arg1;
     246                 :         16 :                 mp_float_t exp_x3 = MICROPY_FLOAT_C_FUN(exp)(x3);
     247                 :         16 :                 lhs_real = exp_x3 * MICROPY_FLOAT_C_FUN(cos)(y3);
     248                 :         16 :                 lhs_imag = exp_x3 * MICROPY_FLOAT_C_FUN(sin)(y3);
     249                 :            :             }
     250                 :            :             break;
     251                 :            :         }
     252                 :            : 
     253                 :       1200 :         case MP_BINARY_OP_EQUAL:
     254   [ +  +  +  + ]:       1200 :             return mp_obj_new_bool(lhs_real == rhs_real && lhs_imag == rhs_imag);
     255                 :            : 
     256                 :            :         default:
     257                 :            :             return MP_OBJ_NULL; // op not supported
     258                 :            :     }
     259                 :        152 :     return mp_obj_new_complex(lhs_real, lhs_imag);
     260                 :            : }
     261                 :            : 
     262                 :            : #endif

Generated by: LCOV version 1.15-5-g462f71d