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