set_d.c
上传用户:qaz666999
上传日期:2022-08-06
资源大小:2570k
文件大小:4k
源码类别:

数学计算

开发平台:

Unix_Linux

  1. /* mpq_set_d(mpq_t q, double d) -- Set q to d without rounding.
  2. Copyright 2000, 2002, 2003 Free Software Foundation, Inc.
  3. This file is part of the GNU MP Library.
  4. The GNU MP Library is free software; you can redistribute it and/or modify
  5. it under the terms of the GNU Lesser General Public License as published by
  6. the Free Software Foundation; either version 3 of the License, or (at your
  7. option) any later version.
  8. The GNU MP Library is distributed in the hope that it will be useful, but
  9. WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
  10. or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
  11. License for more details.
  12. You should have received a copy of the GNU Lesser General Public License
  13. along with the GNU MP Library.  If not, see http://www.gnu.org/licenses/.  */
  14. #include "config.h"
  15. #if HAVE_FLOAT_H
  16. #include <float.h>  /* for DBL_MAX */
  17. #endif
  18. #include "gmp.h"
  19. #include "gmp-impl.h"
  20. #include "longlong.h"
  21. #if LIMBS_PER_DOUBLE > 4
  22.   choke me
  23. #endif
  24. void
  25. mpq_set_d (mpq_ptr dest, double d)
  26. {
  27.   int negative;
  28.   mp_exp_t exp;
  29.   mp_limb_t tp[LIMBS_PER_DOUBLE];
  30.   mp_ptr np, dp;
  31.   mp_size_t nn, dn;
  32.   int c;
  33.   DOUBLE_NAN_INF_ACTION (d,
  34.                          __gmp_invalid_operation (),
  35.                          __gmp_invalid_operation ());
  36.   negative = d < 0;
  37.   d = ABS (d);
  38.   exp = __gmp_extract_double (tp, d);
  39.   /* There are two main version of the conversion.  The `then' arm handles
  40.      numbers with a fractional part, while the `else' arm handles integers.  */
  41. #if LIMBS_PER_DOUBLE == 4
  42.   if (exp <= 1 || (exp == 2 && (tp[0] | tp[1]) != 0))
  43. #endif
  44. #if LIMBS_PER_DOUBLE == 3
  45.   if (exp <= 1 || (exp == 2 && tp[0] != 0))
  46. #endif
  47. #if LIMBS_PER_DOUBLE == 2
  48.   if (exp <= 1)
  49. #endif
  50.     {
  51.       if (d == 0.0)
  52. {
  53.   SIZ(&(dest->_mp_num)) = 0;
  54.   SIZ(&(dest->_mp_den)) = 1;
  55.   PTR(&(dest->_mp_den))[0] = 1;
  56.   return;
  57. }
  58.       dn = -exp;
  59.       MPZ_REALLOC (&(dest->_mp_num), 3);
  60.       np = PTR(&(dest->_mp_num));
  61. #if LIMBS_PER_DOUBLE == 4
  62.       if ((tp[0] | tp[1] | tp[2]) == 0)
  63. np[0] = tp[3], nn = 1;
  64.       else if ((tp[0] | tp[1]) == 0)
  65. np[1] = tp[3], np[0] = tp[2], nn = 2;
  66.       else if (tp[0] == 0)
  67. np[2] = tp[3], np[1] = tp[2], np[0] = tp[1], nn = 3;
  68.       else
  69. np[3] = tp[3], np[2] = tp[2], np[1] = tp[1], np[0] = tp[0], nn = 4;
  70. #endif
  71. #if LIMBS_PER_DOUBLE == 3
  72.       if ((tp[0] | tp[1]) == 0)
  73. np[0] = tp[2], nn = 1;
  74.       else if (tp[0] == 0)
  75. np[1] = tp[2], np[0] = tp[1], nn = 2;
  76.       else
  77. np[2] = tp[2], np[1] = tp[1], np[0] = tp[0], nn = 3;
  78. #endif
  79. #if LIMBS_PER_DOUBLE == 2
  80.       if (tp[0] == 0)
  81. np[0] = tp[1], nn = 1;
  82.       else
  83. np[1] = tp[1], np[0] = tp[0], nn = 2;
  84. #endif
  85.       dn += nn + 1;
  86.       ASSERT_ALWAYS (dn > 0);
  87.       MPZ_REALLOC (&(dest->_mp_den), dn);
  88.       dp = PTR(&(dest->_mp_den));
  89.       MPN_ZERO (dp, dn - 1);
  90.       dp[dn - 1] = 1;
  91.       count_trailing_zeros (c, np[0] | dp[0]);
  92.       if (c != 0)
  93. {
  94.   mpn_rshift (np, np, nn, c);
  95.   nn -= np[nn - 1] == 0;
  96.   mpn_rshift (dp, dp, dn, c);
  97.   dn -= dp[dn - 1] == 0;
  98. }
  99.       SIZ(&(dest->_mp_den)) = dn;
  100.       SIZ(&(dest->_mp_num)) = negative ? -nn : nn;
  101.     }
  102.   else
  103.     {
  104.       nn = exp;
  105.       MPZ_REALLOC (&(dest->_mp_num), nn);
  106.       np = PTR(&(dest->_mp_num));
  107.       switch (nn)
  108.         {
  109. default:
  110.   MPN_ZERO (np, nn - LIMBS_PER_DOUBLE);
  111.   np += nn - LIMBS_PER_DOUBLE;
  112.   /* fall through */
  113. #if LIMBS_PER_DOUBLE == 2
  114. case 2:
  115.   np[1] = tp[1], np[0] = tp[0];
  116.   break;
  117. #endif
  118. #if LIMBS_PER_DOUBLE == 3
  119. case 3:
  120.   np[2] = tp[2], np[1] = tp[1], np[0] = tp[0];
  121.   break;
  122. case 2:
  123.   np[1] = tp[2], np[0] = tp[1];
  124.   break;
  125. #endif
  126. #if LIMBS_PER_DOUBLE == 4
  127. case 4:
  128.   np[3] = tp[3], np[2] = tp[2], np[1] = tp[1], np[0] = tp[0];
  129.   break;
  130. case 3:
  131.   np[2] = tp[3], np[1] = tp[2], np[0] = tp[1];
  132.   break;
  133. case 2:
  134.   np[1] = tp[3], np[0] = tp[2];
  135.   break;
  136. #endif
  137. }
  138.       dp = PTR(&(dest->_mp_den));
  139.       dp[0] = 1;
  140.       SIZ(&(dest->_mp_den)) = 1;
  141.       SIZ(&(dest->_mp_num)) = negative ? -nn : nn;
  142.     }
  143. }