1 /*
2  * Copyright (C) 1987-2008 Sun Microsystems, Inc. All Rights Reserved.
3  * Copyright (C) 2008-2011 Robert Ancell.
4  *
5  * This program is free software: you can redistribute it and/or modify it under
6  * the terms of the GNU General Public License as published by the Free Software
7  * Foundation, either version 2 of the License, or (at your option) any later
8  * version. See http://www.gnu.org/copyleft/gpl.html the full text of the
9  * license.
10  */
11 
12 /*  This maths library is based on the MP multi-precision floating-point
13  *  arithmetic package originally written in FORTRAN by Richard Brent,
14  *  Computer Centre, Australian National University in the 1970's.
15  *
16  *  It has been converted from FORTRAN into C using the freely available
17  *  f2c translator, available via netlib on research.att.com.
18  *
19  *  The subsequently converted C code has then been tidied up, mainly to
20  *  remove any dependencies on the libI77 and libF77 support libraries.
21  *
22  *  FOR A GENERAL DESCRIPTION OF THE PHILOSOPHY AND DESIGN OF MP,
23  *  SEE - R. P. BRENT, A FORTRAN MULTIPLE-PRECISION ARITHMETIC
24  *  PACKAGE, ACM TRANS. MATH. SOFTWARE 4 (MARCH 1978), 57-70.
25  *  SOME ADDITIONAL DETAILS ARE GIVEN IN THE SAME ISSUE, 71-81.
26  *  FOR DETAILS OF THE IMPLEMENTATION, CALLING SEQUENCES ETC. SEE
27  *  THE MP USERS GUIDE.
28  */
29 
30 #ifndef MP_H
31 #define MP_H
32 
33 #include <stdbool.h>
34 #include <stdint.h>
35 #include <glib.h>
36 #include <glib/gi18n.h>
37 #include <mpfr.h>
38 #include <mpc.h>
39 
40 /* If we're not using GNU C, elide __attribute__ */
41 #ifndef __GNUC__
42 #  define  __attribute__(x)  /*NOTHING*/
43 #endif
44 
45 /* Precision of mpfr_t and mpc_t type objects */
46 #define PRECISION 1000
47 
48 typedef struct
49 {
50     mpc_t num;
51 } MPNumber;
52 
53 typedef enum
54 {
55     MP_RADIANS,
56     MP_DEGREES,
57     MP_GRADIANS
58 } MPAngleUnit;
59 
60 /* Returns error string or NULL if no error */
61 // FIXME: Global variable
62 const char  *mp_get_error(void);
63 
64 /* Clear any current error */
65 void        mp_clear_error(void);
66 
67 void        mperr(const char *format, ...) __attribute__((format(printf, 1, 2)));
68 
69 /* Returns initialized MPNumber object */
70 MPNumber    mp_new(void);
71 
72 MPNumber    mp_new_from_unsigned_integer(unsigned long x);
73 
74 MPNumber*   mp_new_ptr(void);
75 
76 void        mp_clear(MPNumber *z);
77 
78 void        mp_free(MPNumber *z);
79 
80 /* Returns:
81  *  0 if x == y
82  * <0 if x < y
83  * >0 if x > y
84  */
85 int    mp_compare(const MPNumber *x, const MPNumber *y);
86 
87 /* Return true if the value is x == 0 */
88 bool   mp_is_zero(const MPNumber *x);
89 
90 /* Return true if x < 0 */
91 bool   mp_is_negative(const MPNumber *x);
92 
93 /* Return true if x is integer */
94 bool   mp_is_integer(const MPNumber *x);
95 
96 /* Return true if x is a positive integer */
97 bool   mp_is_positive_integer(const MPNumber *x);
98 
99 /* Return true if x is a natural number (an integer ≥ 0) */
100 bool   mp_is_natural(const MPNumber *x);
101 
102 /* Return true if x has an imaginary component */
103 bool   mp_is_complex(const MPNumber *x);
104 
105 /* Return true if x == y */
106 bool   mp_is_equal(const MPNumber *x, const MPNumber *y);
107 
108 /* Return true if x ≥ y */
109 bool   mp_is_greater_equal(const MPNumber *x, const MPNumber *y);
110 
111 /* Return true if x > y */
112 bool   mp_is_greater_than(const MPNumber *x, const MPNumber *y);
113 
114 /* Return true if x ≤ y */
115 bool   mp_is_less_equal(const MPNumber *x, const MPNumber *y);
116 
117 /* Return true if x < y */
118 bool   mp_is_less_than(const MPNumber *x, const MPNumber *y);
119 
120 /* Sets z = |x| */
121 void   mp_abs(const MPNumber *x, MPNumber *z);
122 
123 /* Sets z = Arg(x) */
124 void   mp_arg(const MPNumber *x, MPAngleUnit unit, MPNumber *z);
125 
126 /* Sets z = ‾̅x */
127 void   mp_conjugate(const MPNumber *x, MPNumber *z);
128 
129 /* Sets z = Re(x) */
130 void   mp_real_component(const MPNumber *x, MPNumber *z);
131 
132 /* Sets z = Im(x) */
133 void   mp_imaginary_component(const MPNumber *x, MPNumber *z);
134 
135 /* Sets z = −x */
136 void   mp_invert_sign(const MPNumber *x, MPNumber *z);
137 
138 /* Sets z = x + y */
139 void   mp_add(const MPNumber *x, const MPNumber *y, MPNumber *z);
140 
141 /* Sets z = x + y */
142 void   mp_add_integer(const MPNumber *x, long y, MPNumber *z);
143 
144 /* Sets z = x − y */
145 void   mp_subtract(const MPNumber *x, const MPNumber *y, MPNumber *z);
146 
147 /* Sets z = x × y */
148 void   mp_multiply(const MPNumber *x, const MPNumber *y, MPNumber *z);
149 
150 /* Sets z = x × y */
151 void   mp_multiply_integer(const MPNumber *x, long y, MPNumber *z);
152 
153 /* Sets z = x ÷ y */
154 void   mp_divide(const MPNumber *x, const MPNumber *y, MPNumber *z);
155 
156 /* Sets z = x ÷ y */
157 void   mp_divide_integer(const MPNumber *x, long y, MPNumber *z);
158 
159 /* Sets z = 1 ÷ x */
160 void   mp_reciprocal(const MPNumber *, MPNumber *);
161 
162 /* Sets z = sgn(x) */
163 void   mp_sgn(const MPNumber *x, MPNumber *z);
164 
165 void   mp_integer_component(const MPNumber *x, MPNumber *z);
166 
167 /* Sets z = x mod 1 */
168 void   mp_fractional_component(const MPNumber *x, MPNumber *z);
169 
170 /* Sets z = {x} */
171 void   mp_fractional_part(const MPNumber *x, MPNumber *z);
172 
173 /* Sets z = ⌊x⌋ */
174 void   mp_floor(const MPNumber *x, MPNumber *z);
175 
176 /* Sets z = ⌈x⌉ */
177 void   mp_ceiling(const MPNumber *x, MPNumber *z);
178 
179 /* Sets z = [x] */
180 void   mp_round(const MPNumber *x, MPNumber *z);
181 
182 /* Sets z = ln x */
183 void   mp_ln(const MPNumber *x, MPNumber *z);
184 
185 /* Sets z = log_n x */
186 void   mp_logarithm(long n, const MPNumber *x, MPNumber *z);
187 
188 /* Sets z = π */
189 void   mp_get_pi(MPNumber *z);
190 
191 /* Sets z = e */
192 void   mp_get_eulers(MPNumber *z);
193 
194 /* Sets z = i (√−1) */
195 void   mp_get_i(MPNumber *z);
196 
197 /* Sets z = n√x */
198 void   mp_root(const MPNumber *x, long n, MPNumber *z);
199 
200 /* Sets z = √x */
201 void   mp_sqrt(const MPNumber *x, MPNumber *z);
202 
203 /* Sets z = x! */
204 void   mp_factorial(const MPNumber *x, MPNumber *z);
205 
206 /* Sets z = x mod y */
207 void   mp_modulus_divide(const MPNumber *x, const MPNumber *y, MPNumber *z);
208 
209 /* Sets z = x ^ y mod p */
210 void mp_modular_exponentiation(const MPNumber *x, const MPNumber *y, const MPNumber *p, MPNumber *z);
211 
212 /* Sets z = x^y */
213 void   mp_xpowy(const MPNumber *x, const MPNumber *y, MPNumber *z);
214 
215 /* Sets z = x^y */
216 void   mp_xpowy_integer(const MPNumber *x, long y, MPNumber *z);
217 
218 /* Sets z = e^x */
219 void   mp_epowy(const MPNumber *x, MPNumber *z);
220 
221 /* Returns a list of all prime factors in x as MPNumbers */
222 GList* mp_factorize(const MPNumber *x);
223 
224 GList* mp_factorize_unit64 (uint64_t n);
225 
226 /* Sets z = x */
227 void   mp_set_from_mp(const MPNumber *x, MPNumber *z);
228 
229 /* Sets z = x */
230 void   mp_set_from_float(float x, MPNumber *z);
231 
232 /* Sets z = x */
233 void   mp_set_from_double(double x, MPNumber *z);
234 
235 /* Sets z = x */
236 void   mp_set_from_integer(long x, MPNumber *z);
237 
238 /* Sets z = x */
239 void   mp_set_from_unsigned_integer(unsigned long x, MPNumber *z);
240 
241 /* Sets z = numerator ÷ denominator */
242 void   mp_set_from_fraction(long numerator, long denominator, MPNumber *z);
243 
244 /* Sets z = r(cos theta + i sin theta) */
245 void   mp_set_from_polar(const MPNumber *r, MPAngleUnit unit, const MPNumber *theta, MPNumber *z);
246 
247 /* Sets z = x + iy */
248 void   mp_set_from_complex(const MPNumber *x, const MPNumber *y, MPNumber *z);
249 
250 /* Sets z to be a uniform random number in the range [0, 1] */
251 void   mp_set_from_random(MPNumber *z);
252 
253 /* Sets z from a string representation in 'text'.
254  * Returns true on success.
255  */
256 bool   mp_set_from_string(const char *text, int default_base, MPNumber *z);
257 
258 /* Returns x as a native single-precision floating point number */
259 float  mp_to_float(const MPNumber *x);
260 
261 /* Returns x as a native double-precision floating point number */
262 double mp_to_double(const MPNumber *x);
263 
264 /* Returns x as a native integer */
265 long mp_to_integer(const MPNumber *x);
266 
267 /* Returns x as a native unsigned integer */
268 unsigned long mp_to_unsigned_integer(const MPNumber *x);
269 
270 /* Sets z = sin x */
271 void   mp_sin(const MPNumber *x, MPAngleUnit unit, MPNumber *z);
272 
273 /* Sets z = cos x */
274 void   mp_cos(const MPNumber *x, MPAngleUnit unit, MPNumber *z);
275 
276 /* Sets z = tan x */
277 void   mp_tan(const MPNumber *x, MPAngleUnit unit, MPNumber *z);
278 
279 /* Sets z = sin⁻¹ x */
280 void   mp_asin(const MPNumber *x, MPAngleUnit unit, MPNumber *z);
281 
282 /* Sets z = cos⁻¹ x */
283 void   mp_acos(const MPNumber *x, MPAngleUnit unit, MPNumber *z);
284 
285 /* Sets z = tan⁻¹ x */
286 void   mp_atan(const MPNumber *x, MPAngleUnit unit, MPNumber *z);
287 
288 /* Sets z = sinh x */
289 void   mp_sinh(const MPNumber *x, MPNumber *z);
290 
291 /* Sets z = cosh x */
292 void   mp_cosh(const MPNumber *x, MPNumber *z);
293 
294 /* Sets z = tanh x */
295 void   mp_tanh(const MPNumber *x, MPNumber *z);
296 
297 /* Sets z = sinh⁻¹ x */
298 void   mp_asinh(const MPNumber *x, MPNumber *z);
299 
300 /* Sets z = cosh⁻¹ x */
301 void   mp_acosh(const MPNumber *x, MPNumber *z);
302 
303 /* Sets z = tanh⁻¹ x */
304 void   mp_atanh(const MPNumber *x, MPNumber *z);
305 
306 /* Sets z to the value of the error function of x */
307 void   mp_erf(const MPNumber *x, MPNumber *z);
308 
309 /* Sets z to the value of the Riemann Zeta function of x */
310 void   mp_zeta(const MPNumber *x, MPNumber *z);
311 
312 /* Returns true if x is cannot be represented in a binary word of length 'wordlen' */
313 bool   mp_is_overflow(const MPNumber *x, int wordlen);
314 
315 /* Sets z = boolean AND for each bit in x and z */
316 void   mp_and(const MPNumber *x, const MPNumber *y, MPNumber *z);
317 
318 /* Sets z = boolean OR for each bit in x and z */
319 void   mp_or(const MPNumber *x, const MPNumber *y, MPNumber *z);
320 
321 /* Sets z = boolean XOR for each bit in x and z */
322 void   mp_xor(const MPNumber *x, const MPNumber *y, MPNumber *z);
323 
324 /* Sets z = boolean XNOR for each bit in x and z for word of length 'wordlen' */
325 void   mp_xnor(const MPNumber *x, const MPNumber *y, int wordlen, MPNumber *z);
326 
327 /* Sets z = boolean NOT for each bit in x and z for word of length 'wordlen' */
328 void   mp_not(const MPNumber *x, int wordlen, MPNumber *z);
329 
330 /* Sets z = x shifted by 'count' bits.  Positive shift increases the value, negative decreases */
331 void   mp_shift(const MPNumber *x, int count, MPNumber *z);
332 
333 /* Sets z to be the ones complement of x for word of length 'wordlen' */
334 void   mp_ones_complement(const MPNumber *x, int wordlen, MPNumber *z);
335 
336 /* Sets z to be the twos complement of x for word of length 'wordlen' */
337 void   mp_twos_complement(const MPNumber *x, int wordlen, MPNumber *z);
338 
339 void convert_to_radians(const MPNumber *x, MPAngleUnit unit, MPNumber *z);
340 
341 void convert_from_radians(const MPNumber *x, MPAngleUnit unit, MPNumber *z);
342 
343 #endif /* MP_H */
344