1 /* Test file for the beta function
2
3 Copyright 2017-2023 Free Software Foundation, Inc.
4 Contributed by ChemicalDevelopment.
5
6 This file is part of the GNU MPFR Library.
7
8 The GNU MPFR Library is free software; you can redistribute it and/or modify
9 it under the terms of the GNU Lesser General Public License as published by
10 the Free Software Foundation; either version 3 of the License, or (at your
11 option) any later version.
12
13 The GNU MPFR Library is distributed in the hope that it will be useful, but
14 WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
15 or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
16 License for more details.
17
18 You should have received a copy of the GNU Lesser General Public License
19 along with the GNU MPFR Library; see the file COPYING.LESSER. If not, see
20 https://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc.,
21 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA. */
22
23 #include "mpfr-test.h"
24
25 /* TODO: Test the ternary value and the flags. Add tgeneric tests. */
26
27 #define FAILED(p, r, z, w, expected, rnd_mode) do { \
28 printf ("prec=%d, rnd=%s case failed for:", \
29 (int) p, mpfr_print_rnd_mode (rnd_mode)); \
30 printf("\n z ="); \
31 mpfr_out_str (stdout, 2, 0, z, MPFR_RNDN); \
32 printf("\n w ="); \
33 mpfr_out_str (stdout, 2, 0, w, MPFR_RNDN); \
34 printf("\n expected "); \
35 mpfr_out_str (stdout, 2, 0, expected, MPFR_RNDN); \
36 printf("\n got "); \
37 mpfr_out_str (stdout, 2, 0, r, MPFR_RNDN); \
38 printf("\n\n"); \
39 exit (1); \
40 } while (0)
41
42 #define TEST(p, r, z, w, expected) TESTRND(p, r, z, w, expected, MPFR_RNDN)
43
44 #define TESTRND(p, r, z, w, expected, rnd_mode) do { \
45 mpfr_beta (r, z, w, rnd_mode); \
46 if (! SAME_VAL (r, expected)) \
47 FAILED(p, r, z, w, expected, rnd_mode); \
48 } while (0)
49
50 static void
51 test_beta_special (mpfr_prec_t prec)
52 {
53 mpfr_t z, w, r, expect;
54
55 mpfr_init2 (r, prec);
56 mpfr_init2 (z, prec);
57 mpfr_init2 (w, prec);
58 mpfr_init2 (expect, prec);
59
60 mpfr_set_inf (z, 1);
61 mpfr_set_inf (w, 1);
62 mpfr_set_zero (expect, 1);
63 TEST(prec, r, z, w, expect);
64
65 mpfr_set_inf (z, 1);
66 mpfr_set_zero (w, 1);
67 mpfr_set_nan (expect);
68 TEST(prec, r, z, w, expect);
69
70 mpfr_set_inf (z, 1);
71 mpfr_set_zero (w, -1);
72 mpfr_set_nan (expect);
73 TEST(prec, r, z, w, expect);
74
75 mpfr_set_inf (z, 1);
76 mpfr_set_str (w, "-.1e0", 2, MPFR_RNDN);
77 mpfr_set_inf (expect, 1);
78 TEST(prec, r, z, w, expect);
79
80 mpfr_set_inf (z, 1);
81 mpfr_set_str (w, "-1.1e0", 2, MPFR_RNDN);
82 mpfr_set_inf (expect, -1);
83 TEST(prec, r, z, w, expect);
84
85 mpfr_set_inf (z, 1);
86 mpfr_set_str (w, "-1e0", 2, MPFR_RNDN);
87 mpfr_set_nan (expect);
88 TEST(prec, r, z, w, expect);
89
90 mpfr_set_inf (z, 1);
91 mpfr_set_str (w, "-2e0", 2, MPFR_RNDN);
92 mpfr_set_nan (expect);
93 TEST(prec, r, z, w, expect);
94
95 if (prec > 81)
96 {
97 mpfr_set_inf (z, 1);
98 mpfr_set_str (w, "-1e80", 2, MPFR_RNDN);
99 mpfr_set_nan (expect);
100 TEST(prec, r, z, w, expect);
101
102 mpfr_set_inf (z, 1);
103 mpfr_set_str (w, "-1e80", 2, MPFR_RNDN);
104 mpfr_sub_d (w, w, .1, MPFR_RNDN);
105 mpfr_set_inf (expect, 1);
106 TEST(prec, r, z, w, expect);
107
108 mpfr_set_str (w, "-1e80", 2, MPFR_RNDN);
109 mpfr_sub_d (w, w, 1.1, MPFR_RNDN);
110 mpfr_set_inf (expect, -1);
111 TEST(prec, r, z, w, expect);
112 }
113
114 mpfr_set_str (z, "1.1e0", 2, MPFR_RNDN);
115 mpfr_set_inf (w, -1);
116 mpfr_set_nan (expect);
117 TEST(prec, r, z, w, expect);
118
119 mpfr_set_str (z, "11e0", 2, MPFR_RNDN);
120 mpfr_set_inf (w, -1);
121 mpfr_set_zero (expect, -1);
122 TEST(prec, r, z, w, expect);
123
124 mpfr_set_str (z, "10e0", 2, MPFR_RNDN);
125 mpfr_set_inf (w, -1);
126 mpfr_set_zero (expect, 1);
127 TEST(prec, r, z, w, expect);
128
129 mpfr_set_inf (z, -1);
130 mpfr_set_inf (w, -1);
131 mpfr_set_nan (expect);
132 TEST(prec, r, z, w, expect);
133
134 mpfr_clear (r);
135 mpfr_clear (z);
136 mpfr_clear (w);
137 mpfr_clear (expect);
138 }
139
140 static void
141 test_beta_2exp (mpfr_prec_t prec, int trials, int spread)
142 {
143 mpfr_t r, z, w, expect;
144 int i;
145
146 mpfr_init2 (r, prec);
147 mpfr_init2 (z, prec);
148 mpfr_init2 (w, prec);
149 mpfr_init2 (expect, prec);
150 for (i = -(spread*trials)/2; spread*i < trials / 2; i += spread)
151 {
152 mpfr_set_si_2exp (z, 1, i, MPFR_RNDN);
153 mpfr_set_ui (w, 1, MPFR_RNDN);
154 mpfr_set_si_2exp (expect, 1, -i, MPFR_RNDN);
155
156 TEST(prec, r, z, w, expect);
157 }
158
159 mpfr_clear (r);
160 mpfr_clear (z);
161 mpfr_clear (w);
162 mpfr_clear (expect);
163 }
164
165 /*
166 Tests values such that z and w are not integers, but (z+w) is.
167
168 An example that was given:
169 beta(-.3, -1.7) = gamma(-0.3)*gamma(-1.7)/gamma(-2)
170
171 Sage gives this as 0, and Lefevre said that we should return +0
172
173 */
174 static void
175 test_beta_zw_sum_int (mpfr_prec_t prec)
176 {
177 mpfr_t r, z, w, expect;
178 int sum;
179
180 if (prec < 4)
181 prec = 4;
182
183 mpfr_init2 (r, prec);
184 mpfr_init2 (z, prec);
185 mpfr_init2 (w, prec);
186 mpfr_init2 (expect, prec);
187
188 sum = -3;
189 mpfr_set_str (z, "-1.1e0", 2, MPFR_RNDN);
190 mpfr_si_sub (w, sum, z, MPFR_RNDN);
191 mpfr_set_zero (expect, 1);
192 TEST(prec, r, z, w, expect);
193
194 sum = -12;
195 mpfr_set_str (z, "-1.101e0", 2, MPFR_RNDN);
196 mpfr_si_sub (w, sum, z, MPFR_RNDN);
197 mpfr_set_zero (expect, 1);
198 TEST(prec, r, z, w, expect);
199
200 sum = -1;
201 mpfr_set_str (z, "-.11e0", 2, MPFR_RNDN);
202 mpfr_si_sub (w, sum, z, MPFR_RNDN);
203 mpfr_set_zero (expect, 1);
204 TEST(prec, r, z, w, expect);
205
206 sum = -13;
207 mpfr_set_str (z, "-.001e0", 2, MPFR_RNDN);
208 mpfr_si_sub (w, sum, z, MPFR_RNDN);
209 mpfr_set_zero (expect, 1);
210 TEST(prec, r, z, w, expect);
211
212 mpfr_clear (r);
213 mpfr_clear (z);
214 mpfr_clear (w);
215 mpfr_clear (expect);
216 }
217
218
219 static void
220 test_beta_hardcoded (mpfr_prec_t prec)
221 {
222 mpfr_t r, z, w, expect;
223 mpfr_prec_t oprec = 1;
224
225 if (prec < 10)
226 prec = 10;
227
228 mpfr_init2 (z, prec);
229 mpfr_init2 (w, prec);
230 mpfr_init2 (r, oprec);
231 mpfr_init2 (expect, oprec);
232
233 mpfr_set_ui (z, 3, MPFR_RNDN);
234 mpfr_set_ui (w, 3, MPFR_RNDN);
235 mpfr_set_str (expect, "1e-5", 2, MPFR_RNDN);
236 TESTRND(prec, r, z, w, expect, MPFR_RNDN);
237
238 mpfr_set_str (expect, "1.1e-5", 2, MPFR_RNDN);
239 TESTRND(prec, r, z, w, expect, MPFR_RNDU);
240
241 mpfr_set_str (expect, "1e-5", 2, MPFR_RNDN);
242 TESTRND(prec, r, z, w, expect, MPFR_RNDD);
243
244 mpfr_set_ui (z, 5, MPFR_RNDN);
245 mpfr_set_ui (w, 27, MPFR_RNDN);
246 mpfr_set_str (expect, "1e-20", 2, MPFR_RNDN);
247 TESTRND(prec, r, z, w, expect, MPFR_RNDN);
248
249 mpfr_set_str (expect, "1e-19", 2, MPFR_RNDN);
250 TESTRND(prec, r, z, w, expect, MPFR_RNDU);
251
252 mpfr_set_ui (z, 5, MPFR_RNDN);
253 mpfr_set_ui (w, 27, MPFR_RNDN);
254 mpfr_set_str (expect, "1e-20", 2, MPFR_RNDN);
255 TESTRND(prec, r, z, w, expect, MPFR_RNDN);
256
257
258 mpfr_set_ui (z, 121, MPFR_RNDN);
259 mpfr_set_ui (w, 2, MPFR_RNDN);
260 mpfr_set_str (expect, "1e-14", 2, MPFR_RNDN);
261 TESTRND(prec, r, z, w, expect, MPFR_RNDN);
262
263 mpfr_set_ui (z, 121, MPFR_RNDN);
264 mpfr_set_ui (w, 151, MPFR_RNDN);
265 mpfr_set_str (expect, "1e-271", 2, MPFR_RNDN);
266 TESTRND(prec, r, z, w, expect, MPFR_RNDN);
267
268 mpfr_set_str (expect, "1e-272", 2, MPFR_RNDN);
269 TESTRND(prec, r, z, w, expect, MPFR_RNDD);
270
271 mpfr_set_str (expect, "1e-271", 2, MPFR_RNDN);
272 TESTRND(prec, r, z, w, expect, MPFR_RNDU);
273
274 mpfr_clear (r);
275 mpfr_clear (z);
276 mpfr_clear (w);
277 mpfr_clear (expect);
278 }
279
280 /* makes sure beta(a, b) = beta(b, a) */
281 static void
282 test_beta_refl (mpfr_prec_t prec, mpfr_rnd_t rnd_mode)
283 {
284 mpfr_t r, z, w, expect;
285
286 mpfr_init2 (z, prec);
287 mpfr_init2 (w, prec);
288 mpfr_init2 (r, prec);
289 mpfr_init2 (expect, prec);
290
291 mpfr_set_ui (z, 3, MPFR_RNDN);
292 mpfr_set_ui (w, 3, MPFR_RNDN);
293 mpfr_beta (expect, w, z, rnd_mode);
294 TESTRND(prec, r, z, w, expect, rnd_mode);
295
296 mpfr_set_ui (z, 5, MPFR_RNDN);
297 mpfr_set_ui (w, 100, MPFR_RNDN);
298 mpfr_beta (expect, w, z, rnd_mode);
299 TESTRND(prec, r, z, w, expect, rnd_mode);
300
301 mpfr_set_nan (z);
302 mpfr_set_ui (w, 100, MPFR_RNDN);
303 mpfr_beta (expect, w, z, rnd_mode);
304 TESTRND(prec, r, z, w, expect, rnd_mode);
305
306 mpfr_set_nan (z);
307 mpfr_set_ui (w, 1, MPFR_RNDN);
308 mpfr_beta (expect, w, z, rnd_mode);
309 TESTRND(prec, r, z, w, expect, rnd_mode);
310
311 mpfr_set_nan (z);
312 mpfr_set_nan (w);
313 mpfr_beta (expect, w, z, rnd_mode);
314 TESTRND(prec, r, z, w, expect, rnd_mode);
315 mpfr_set_nan (z);
316 mpfr_set_nan (w);
317 mpfr_beta (expect, w, z, rnd_mode);
318 TESTRND(prec, r, z, w, expect, rnd_mode);
319
320 mpfr_clear (r);
321 mpfr_clear (z);
322 mpfr_clear (w);
323 mpfr_clear (expect);
324 }
325
326 #define TEST_FUNCTION mpfr_beta
327 #define TWO_ARGS
328 #define TEST_RANDOM_EMIN -16
329 #define TEST_RANDOM_EMAX 16
330 #include "tgeneric.c"
331
332 int
333 main (void)
334 {
335 tests_start_mpfr ();
336
337 test_beta_special (10);
338 test_beta_special (100);
339 test_beta_special (1000);
340
341 test_beta_2exp (1, 10, 1);
342 test_beta_2exp (100, 40, 3);
343
344 test_beta_hardcoded (10);
345 test_beta_hardcoded (100);
346
347 test_beta_refl (1, MPFR_RNDN);
348 test_beta_refl (100, MPFR_RNDD);
349
350 test_beta_zw_sum_int (10);
351 test_beta_zw_sum_int (100);
352
353 test_generic (MPFR_PREC_MIN, 100, 20);
354
355 tests_end_mpfr ();
356 return 0;
357 }