/contrib/mpc/src/get_x.c

https://bitbucket.org/cooljeanius/dragonflybsd · C · 236 lines · 152 code · 47 blank · 37 comment · 22 complexity · a1c844b2298a075f0821ebf3099a0f3e MD5 · raw file

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