Home | History | Annotate | Line # | Download | only in src
      1      1.1  mrg /* mpc_get_dc, mpc_get_ldc -- Transform mpc number into C complex number
      2      1.1  mrg    mpc_get_str -- Convert a complex number into a string.
      3      1.1  mrg 
      4  1.1.1.3  mrg Copyright (C) 2009, 2010, 2011, 2020, 2022 INRIA
      5      1.1  mrg 
      6      1.1  mrg This file is part of GNU MPC.
      7      1.1  mrg 
      8      1.1  mrg GNU MPC is free software; you can redistribute it and/or modify it under
      9      1.1  mrg the terms of the GNU Lesser General Public License as published by the
     10      1.1  mrg Free Software Foundation; either version 3 of the License, or (at your
     11      1.1  mrg option) any later version.
     12      1.1  mrg 
     13      1.1  mrg GNU MPC is distributed in the hope that it will be useful, but WITHOUT ANY
     14      1.1  mrg WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
     15      1.1  mrg FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
     16      1.1  mrg more details.
     17      1.1  mrg 
     18      1.1  mrg You should have received a copy of the GNU Lesser General Public License
     19      1.1  mrg along with this program. If not, see http://www.gnu.org/licenses/ .
     20      1.1  mrg */
     21      1.1  mrg 
     22      1.1  mrg #include "config.h"
     23      1.1  mrg 
     24      1.1  mrg #ifdef HAVE_COMPLEX_H
     25      1.1  mrg #include <complex.h>
     26      1.1  mrg #endif
     27      1.1  mrg 
     28      1.1  mrg #include <locale.h>
     29      1.1  mrg #include <stdio.h> /* for sprintf, fprintf */
     30      1.1  mrg #include <ctype.h>
     31      1.1  mrg #include <string.h>
     32      1.1  mrg #include "mpc-impl.h"
     33      1.1  mrg 
     34      1.1  mrg #ifdef HAVE_COMPLEX_H
     35      1.1  mrg double _Complex
     36      1.1  mrg mpc_get_dc (mpc_srcptr op, mpc_rnd_t rnd) {
     37      1.1  mrg    return I * mpfr_get_d (mpc_imagref (op), MPC_RND_IM (rnd))
     38      1.1  mrg           + mpfr_get_d (mpc_realref (op), MPC_RND_RE (rnd));
     39      1.1  mrg }
     40      1.1  mrg 
     41      1.1  mrg long double _Complex
     42      1.1  mrg mpc_get_ldc (mpc_srcptr op, mpc_rnd_t rnd) {
     43      1.1  mrg    return I * mpfr_get_ld (mpc_imagref (op), MPC_RND_IM (rnd))
     44      1.1  mrg           + mpfr_get_ld (mpc_realref (op), MPC_RND_RE (rnd));
     45      1.1  mrg }
     46      1.1  mrg #endif
     47      1.1  mrg 
     48      1.1  mrg 
     49      1.1  mrg /* Code for mpc_get_str. The output format is "(real imag)", the decimal point
     50      1.1  mrg    of the locale is used. */
     51      1.1  mrg 
     52      1.1  mrg /* mpfr_prec_t can be either int or long int */
     53  1.1.1.2  mrg #if (_MPFR_EXP_FORMAT == 2)
     54      1.1  mrg #define MPC_EXP_FORMAT_SPEC "i"
     55  1.1.1.2  mrg #elif (_MPFR_EXP_FORMAT == 3)
     56      1.1  mrg #define MPC_EXP_FORMAT_SPEC "li"
     57      1.1  mrg #else
     58  1.1.1.2  mrg #error "value of _MPFR_EXP_FORMAT not supported"
     59      1.1  mrg #endif
     60      1.1  mrg 
     61      1.1  mrg static char *
     62      1.1  mrg pretty_zero (mpfr_srcptr zero)
     63      1.1  mrg {
     64      1.1  mrg   char *pretty;
     65      1.1  mrg 
     66      1.1  mrg   pretty = mpc_alloc_str (3);
     67      1.1  mrg 
     68      1.1  mrg   pretty[0] = mpfr_signbit (zero) ? '-' : '+';
     69      1.1  mrg   pretty[1] = '0';
     70      1.1  mrg   pretty[2] = '\0';
     71      1.1  mrg 
     72      1.1  mrg   return pretty;
     73      1.1  mrg }
     74      1.1  mrg 
     75      1.1  mrg static char *
     76  1.1.1.2  mrg prettify (const char *str, const mpfr_exp_t expo, int base, int special)
     77      1.1  mrg {
     78      1.1  mrg   size_t sz;
     79      1.1  mrg   char *pretty;
     80      1.1  mrg   char *p;
     81      1.1  mrg   const char *s;
     82  1.1.1.2  mrg   mpfr_exp_t x;
     83      1.1  mrg   int sign;
     84      1.1  mrg 
     85      1.1  mrg   sz = strlen (str) + 1; /* + terminal '\0' */
     86      1.1  mrg 
     87      1.1  mrg   if (special)
     88      1.1  mrg     {
     89      1.1  mrg       /* special number: nan or inf */
     90      1.1  mrg       pretty = mpc_alloc_str (sz);
     91      1.1  mrg       strcpy (pretty, str);
     92      1.1  mrg 
     93      1.1  mrg       return pretty;
     94      1.1  mrg     }
     95      1.1  mrg 
     96      1.1  mrg   /* regular number */
     97      1.1  mrg 
     98      1.1  mrg   sign = (str[0] == '-' || str[0] == '+');
     99      1.1  mrg 
    100      1.1  mrg   x = expo - 1; /* expo is the exponent value with decimal point BEFORE
    101  1.1.1.2  mrg                    the first digit, we want decimal point AFTER the first
    102      1.1  mrg                    digit */
    103      1.1  mrg   if (base == 16)
    104  1.1.1.2  mrg     x *= 4; /* the output exponent is a binary exponent */
    105      1.1  mrg 
    106      1.1  mrg   ++sz; /* + decimal point */
    107      1.1  mrg 
    108      1.1  mrg   if (x != 0)
    109      1.1  mrg     {
    110      1.1  mrg       /* augment sz with the size needed for an exponent written in base
    111      1.1  mrg          ten */
    112  1.1.1.2  mrg       mpfr_exp_t xx;
    113      1.1  mrg 
    114      1.1  mrg       sz += 3; /* + exponent char + sign + 1 digit */
    115      1.1  mrg 
    116      1.1  mrg       if (x < 0)
    117      1.1  mrg         {
    118      1.1  mrg           /* avoid overflow when changing sign (assuming that, for the
    119  1.1.1.2  mrg              mpfr_exp_t type, (max value) is greater than (- min value / 10)) */
    120      1.1  mrg           if (x < -10)
    121      1.1  mrg             {
    122      1.1  mrg               xx = - (x / 10);
    123      1.1  mrg               sz++;
    124      1.1  mrg             }
    125      1.1  mrg           else
    126      1.1  mrg             xx = -x;
    127      1.1  mrg         }
    128      1.1  mrg       else
    129      1.1  mrg         xx = x;
    130      1.1  mrg 
    131      1.1  mrg       /* compute sz += floor(log(expo)/log(10)) without using libm
    132      1.1  mrg          functions */
    133      1.1  mrg       while (xx > 9)
    134      1.1  mrg         {
    135      1.1  mrg           sz++;
    136      1.1  mrg           xx /= 10;
    137      1.1  mrg         }
    138      1.1  mrg     }
    139      1.1  mrg 
    140      1.1  mrg   pretty = mpc_alloc_str (sz);
    141      1.1  mrg   p = pretty;
    142      1.1  mrg 
    143      1.1  mrg   /* 1. optional sign plus first digit */
    144      1.1  mrg   s = str;
    145      1.1  mrg   *p++ = *s++;
    146      1.1  mrg   if (sign)
    147      1.1  mrg     *p++ = *s++;
    148      1.1  mrg 
    149      1.1  mrg   /* 2. decimal point */
    150      1.1  mrg #ifdef HAVE_LOCALECONV
    151      1.1  mrg   *p++ = *localeconv ()->decimal_point;
    152      1.1  mrg #else
    153      1.1  mrg   *p++ = '.';
    154      1.1  mrg #endif
    155      1.1  mrg   *p = '\0';
    156      1.1  mrg 
    157      1.1  mrg   /* 3. other significant digits */
    158      1.1  mrg   strcat (pretty, s);
    159      1.1  mrg 
    160      1.1  mrg   /* 4. exponent (in base ten) */
    161      1.1  mrg   if (x == 0)
    162      1.1  mrg     return pretty;
    163      1.1  mrg 
    164      1.1  mrg   p = pretty + strlen (str) + 1;
    165      1.1  mrg 
    166      1.1  mrg   switch (base)
    167      1.1  mrg     {
    168      1.1  mrg     case 10:
    169      1.1  mrg       *p++ = 'e';
    170      1.1  mrg       break;
    171      1.1  mrg     case 2:
    172      1.1  mrg     case 16:
    173      1.1  mrg       *p++ = 'p';
    174      1.1  mrg       break;
    175      1.1  mrg     default:
    176      1.1  mrg       *p++ = '@';
    177      1.1  mrg     }
    178      1.1  mrg 
    179      1.1  mrg   *p = '\0';
    180      1.1  mrg 
    181  1.1.1.3  mrg   sprintf (p, "%+" MPC_EXP_FORMAT_SPEC, x);
    182      1.1  mrg 
    183      1.1  mrg   return pretty;
    184      1.1  mrg }
    185      1.1  mrg 
    186      1.1  mrg static char *
    187      1.1  mrg get_pretty_str (const int base, const size_t n, mpfr_srcptr x, mpfr_rnd_t rnd)
    188      1.1  mrg {
    189  1.1.1.2  mrg   mpfr_exp_t expo;
    190      1.1  mrg   char *ugly;
    191      1.1  mrg   char *pretty;
    192      1.1  mrg 
    193      1.1  mrg   if (mpfr_zero_p (x))
    194      1.1  mrg     return pretty_zero (x);
    195      1.1  mrg 
    196      1.1  mrg   ugly = mpfr_get_str (NULL, &expo, base, n, x, rnd);
    197      1.1  mrg   MPC_ASSERT (ugly != NULL);
    198      1.1  mrg   pretty = prettify (ugly, expo, base, !mpfr_number_p (x));
    199      1.1  mrg   mpfr_free_str (ugly);
    200      1.1  mrg 
    201      1.1  mrg   return pretty;
    202      1.1  mrg }
    203      1.1  mrg 
    204      1.1  mrg char *
    205      1.1  mrg mpc_get_str (int base, size_t n, mpc_srcptr op, mpc_rnd_t rnd)
    206      1.1  mrg {
    207      1.1  mrg   size_t needed_size;
    208      1.1  mrg   char *real_str;
    209      1.1  mrg   char *imag_str;
    210      1.1  mrg   char *complex_str = NULL;
    211      1.1  mrg 
    212      1.1  mrg   if (base < 2 || base > 36)
    213      1.1  mrg     return NULL;
    214      1.1  mrg 
    215      1.1  mrg   real_str = get_pretty_str (base, n, mpc_realref (op), MPC_RND_RE (rnd));
    216      1.1  mrg   imag_str = get_pretty_str (base, n, mpc_imagref (op), MPC_RND_IM (rnd));
    217      1.1  mrg 
    218      1.1  mrg   needed_size = strlen (real_str) + strlen (imag_str) + 4;
    219      1.1  mrg 
    220      1.1  mrg   complex_str = mpc_alloc_str (needed_size);
    221      1.1  mrg MPC_ASSERT (complex_str != NULL);
    222      1.1  mrg 
    223      1.1  mrg   strcpy (complex_str, "(");
    224      1.1  mrg   strcat (complex_str, real_str);
    225      1.1  mrg   strcat (complex_str, " ");
    226      1.1  mrg   strcat (complex_str, imag_str);
    227      1.1  mrg   strcat (complex_str, ")");
    228      1.1  mrg 
    229      1.1  mrg   mpc_free_str (real_str);
    230      1.1  mrg   mpc_free_str (imag_str);
    231      1.1  mrg 
    232      1.1  mrg   return complex_str;
    233      1.1  mrg }
    234