xref: /netbsd-src/external/lgpl3/mpfr/dist/tests/tbeta.c (revision ba125506a622fe649968631a56eba5d42ff57863)
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
test_beta_special(mpfr_prec_t prec)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
test_beta_2exp(mpfr_prec_t prec,int trials,int spread)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
test_beta_zw_sum_int(mpfr_prec_t prec)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
test_beta_hardcoded(mpfr_prec_t prec)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
test_beta_refl(mpfr_prec_t prec,mpfr_rnd_t rnd_mode)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
main(void)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 }
358