1 /* mpcheck-float128 -- compare mpc functions against "__float128 complex"
2 from the GNU libc implementation
3
4 Copyright (C) 2020 INRIA
5
6 This file is part of GNU MPC.
7
8 GNU MPC is free software; you can redistribute it and/or modify it under
9 the terms of the GNU Lesser General Public License as published by the
10 Free Software Foundation; either version 3 of the License, or (at your
11 option) any later version.
12
13 GNU MPC is distributed in the hope that it will be useful, but WITHOUT ANY
14 WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
15 FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
16 more details.
17
18 You should have received a copy of the GNU Lesser General Public License
19 along with this program. If not, see http://www.gnu.org/licenses/ .
20 */
21
22 /* the GNU libc provides the following functions (as of 2.31),
23 with 'f' suffix for the float/binary32 version, with no suffix
24 for the double/binary64 version, with 'l' suffix for the long double
25 version, and with 'f128' suffix for the __float128 version:
26
27 cabs casinh cexp csinh
28 cacos catan clog csqrt
29 cacosh catanh clog10 ctan
30 carg ccos cpow ctanh
31 casin ccosh csin
32 */
33
34 #define _GNU_SOURCE /* for clog10 */
35 #include <stdio.h>
36 #include <stdlib.h>
37 #include <string.h>
38 #include <complex.h>
39 #include <sys/types.h>
40 #include <unistd.h>
41 #include <assert.h>
42 #define MPFR_WANT_FLOAT128
43 #include "mpc.h"
44 #ifdef __GNUC__
45 #include <gnu/libc-version.h>
46 #endif
47
48 #define PRECISION 113
49 #define EMAX 16384
50 #define TYPE _Float128
51 #define SUFFIX f128
52
53 #define mpfr_set_type mpfr_set_float128
54
55 static TYPE complex
mpc_get_type(mpc_t z,mpc_rnd_t rnd)56 mpc_get_type (mpc_t z, mpc_rnd_t rnd)
57 {
58 TYPE x, y;
59 /* there is no mpc_get_float128c function */
60 x = mpfr_get_float128 (mpc_realref (z), MPC_RND_RE(rnd));
61 y = mpfr_get_float128 (mpc_imagref (z), MPC_RND_IM(rnd));
62 return x + I * y;
63 }
64
65 static int
mpc_set_type(mpc_t x,TYPE complex y,mpc_rnd_t rnd)66 mpc_set_type (mpc_t x, TYPE complex y, mpc_rnd_t rnd)
67 {
68 /* there is no mpc_set_float128c function */
69 mpfr_set_float128 (mpc_realref (x), crealf128 (y), MPC_RND_RE(rnd));
70 mpfr_set_float128 (mpc_imagref (x), cimagf128 (y), MPC_RND_IM(rnd));
71 }
72
73 gmp_randstate_t state;
74 mpz_t expz; /* global variable used in mpcheck_random */
75 unsigned long seed = 0;
76 int verbose = 0;
77 mpfr_exp_t emin, emax;
78
79 #include "mpcheck-common.c"
80
81 #define FOO add
82 #define CFOO(x,y) (x+y)
83 #include "mpcheck-template3.c"
84
85 #define FOO sub
86 #define CFOO(x,y) (x-y)
87 #include "mpcheck-template3.c"
88
89 #define FOO mul
90 #define CFOO(x,y) (x*y)
91 #include "mpcheck-template3.c"
92
93 #define FOO div
94 #define CFOO(x,y) (x/y)
95 #include "mpcheck-template3.c"
96
97 #define FOO pow
98 #include "mpcheck-template3.c"
99
100 #define FOO abs
101 #include "mpcheck-template2.c"
102
103 #define FOO arg
104 #include "mpcheck-template2.c"
105
106 #define FOO sqrt
107 #include "mpcheck-template1.c"
108
109 #define FOO acos
110 #include "mpcheck-template1.c"
111
112 #define FOO acosh
113 #include "mpcheck-template1.c"
114
115 #define FOO asin
116 #include "mpcheck-template1.c"
117
118 #define FOO asinh
119 #include "mpcheck-template1.c"
120
121 #define FOO atan
122 #include "mpcheck-template1.c"
123
124 #define FOO atanh
125 #include "mpcheck-template1.c"
126
127 #define FOO cos
128 #include "mpcheck-template1.c"
129
130 #define FOO cosh
131 #include "mpcheck-template1.c"
132
133 #define FOO exp
134 #include "mpcheck-template1.c"
135
136 #define FOO log
137 #include "mpcheck-template1.c"
138
139 #define FOO log10
140 #include "mpcheck-template1.c"
141
142 #define FOO sin
143 #include "mpcheck-template1.c"
144
145 #define FOO sinh
146 #include "mpcheck-template1.c"
147
148 #define FOO tan
149 #include "mpcheck-template1.c"
150
151 #define FOO tanh
152 #include "mpcheck-template1.c"
153
154 int
main(int argc,char * argv[])155 main (int argc, char *argv[])
156 {
157 mpfr_prec_t p = PRECISION; /* precision of 'double' */
158 unsigned long n = 1000000; /* default number of random tests per function */
159
160 while (argc >= 2 && argv[1][0] == '-')
161 {
162 if (argc >= 3 && strcmp (argv[1], "-p") == 0)
163 {
164 p = atoi (argv[2]);
165 argc -= 2;
166 argv += 2;
167 }
168 else if (argc >= 3 && strcmp (argv[1], "-seed") == 0)
169 {
170 seed = atoi (argv[2]);
171 argc -= 2;
172 argv += 2;
173 }
174 else if (argc >= 3 && strcmp (argv[1], "-num") == 0)
175 {
176 n = atoi (argv[2]);
177 argc -= 2;
178 argv += 2;
179 }
180 else if (strcmp (argv[1], "-v") == 0)
181 {
182 verbose ++;
183 argc --;
184 argv ++;
185 }
186 else if (strcmp (argv[1], "-check") == 0)
187 {
188 recheck = 1;
189 argc --;
190 argv ++;
191 }
192 else
193 {
194 fprintf (stderr, "Unknown option %s\n", argv[1]);
195 exit (1);
196 }
197 }
198
199 /* set exponent range */
200 emin = -EMAX - 64 + 4; /* should be -16444 like for long double */
201 emax = EMAX;
202 mpfr_set_emin (emin);
203 mpfr_set_emax (emax);
204
205 gmp_randinit_default (state);
206 mpz_init (expz);
207
208 printf ("Using GMP %s, MPFR %s\n", gmp_version, mpfr_get_version ());
209
210 #ifdef __GNUC__
211 printf ("GNU libc version: %s\n", gnu_get_libc_version ());
212 printf ("GNU libc release: %s\n", gnu_get_libc_release ());
213 #endif
214
215 if (seed == 0)
216 seed = getpid ();
217 printf ("Using random seed %lu\n", seed);
218
219 /* (complex,complex) -> complex */
220 test_add (p, n);
221 test_sub (p, n);
222 test_mul (p, n);
223 test_div (p, n);
224 test_pow (p, n);
225
226 /* complex -> real */
227 test_abs (p, n);
228 test_arg (p, n);
229
230 /* complex -> complex */
231 test_sqrt (p, n);
232 test_acos (p, n);
233 test_acosh (p, n);
234 test_asin (p, n);
235 test_asinh (p, n);
236 test_atan (p, n);
237 test_atanh (p, n);
238 test_cos (p, n);
239 test_cosh (p, n);
240 test_exp (p, n);
241 test_log (p, n);
242 test_log10 (p, n);
243 test_sin (p, n);
244 test_sinh (p, n);
245 test_tan (p, n);
246 test_tanh (p, n);
247
248 gmp_randclear (state);
249 mpz_clear (expz);
250
251 report_maximal_errors ();
252
253 return 0;
254 }
255