xref: /netbsd-src/external/lgpl3/mpfr/dist/tests/tlog.c (revision ec6772edaf0cdcb5f52a48f4aca5e33a8fb8ecfd)
1 /* Test file for mpfr_log.
2 
3 Copyright 1999, 2001-2023 Free Software Foundation, Inc.
4 Contributed by the AriC and Caramba projects, INRIA.
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 #ifdef CHECK_EXTERNAL
26 static int
test_log(mpfr_ptr a,mpfr_srcptr b,mpfr_rnd_t rnd_mode)27 test_log (mpfr_ptr a, mpfr_srcptr b, mpfr_rnd_t rnd_mode)
28 {
29   int res;
30   int ok = rnd_mode == MPFR_RNDN && mpfr_number_p (b) && mpfr_get_prec (a)>=53;
31   if (ok)
32     {
33       mpfr_print_raw (b);
34     }
35   res = mpfr_log (a, b, rnd_mode);
36   if (ok)
37     {
38       printf (" ");
39       mpfr_print_raw (a);
40       printf ("\n");
41     }
42   return res;
43 }
44 #else
45 #define test_log mpfr_log
46 #endif
47 
48 static void
check2(const char * as,mpfr_rnd_t rnd_mode,const char * res1s)49 check2 (const char *as, mpfr_rnd_t rnd_mode, const char *res1s)
50 {
51   mpfr_t ta, tres;
52 
53   mpfr_inits2 (53, ta, tres, (mpfr_ptr) 0);
54   mpfr_set_str1 (ta, as);
55   test_log (tres, ta, rnd_mode);
56 
57   if (mpfr_cmp_str1 (tres, res1s))
58     {
59       printf ("mpfr_log failed for    a=%s, rnd_mode=%s\n",
60               as, mpfr_print_rnd_mode (rnd_mode));
61       printf ("correct result is        %s\n mpfr_log gives          ",
62               res1s);
63       mpfr_out_str(stdout, 10, 0, tres, MPFR_RNDN);
64       printf ("\n");
65       exit (1);
66     }
67   mpfr_clears (ta, tres, (mpfr_ptr) 0);
68 }
69 
70 static void
check3(char * s,unsigned long prec,mpfr_rnd_t rnd)71 check3 (char *s, unsigned long prec, mpfr_rnd_t rnd)
72 {
73   mpfr_t x, y;
74 
75   mpfr_init2 (x, prec);
76   mpfr_init2 (y, prec);
77   mpfr_set_str (x, s, 10, rnd);
78   test_log (y, x, rnd);
79   mpfr_out_str (stdout, 10, 0, y, rnd);
80   puts ("");
81   mpfr_clear (x);
82   mpfr_clear (y);
83 }
84 
85 /* examples from Jean-Michel Muller and Vincent Lefevre
86    Cf http://perso.ens-lyon.fr/jean-michel.muller/Intro-to-TMD.htm
87 */
88 
89 static void
check_worst_cases(void)90 check_worst_cases (void)
91 {
92   check2("1.00089971802309629645", MPFR_RNDD, "8.99313519443722736088e-04");
93   check2("1.00089971802309629645", MPFR_RNDN, "8.99313519443722844508e-04");
94   check2("1.00089971802309629645", MPFR_RNDU, "8.99313519443722844508e-04");
95 
96   check2("1.01979300812244555452", MPFR_RNDD, "1.95996734891603630047e-02");
97   check2("1.01979300812244555452", MPFR_RNDN, "1.95996734891603664741e-02");
98   check2("1.01979300812244555452", MPFR_RNDU, "1.95996734891603664741e-02");
99 
100   check2("1.02900871924604464525", MPFR_RNDD, "2.85959303301472726744e-02");
101   check2("1.02900871924604464525", MPFR_RNDN, "2.85959303301472761438e-02");
102   check2("1.02900871924604464525", MPFR_RNDU, "2.85959303301472761438e-02");
103 
104   check2("1.27832870030418943585", MPFR_RNDD, "2.45553521871417795852e-01");
105   check2("1.27832870030418943585", MPFR_RNDN, "2.45553521871417823608e-01");
106   check2("1.27832870030418943585", MPFR_RNDU, "2.45553521871417823608e-01");
107 
108   check2("1.31706530746788241792", MPFR_RNDD, "2.75406009586277422674e-01");
109   check2("1.31706530746788241792", MPFR_RNDN, "2.75406009586277478185e-01");
110   check2("1.31706530746788241792", MPFR_RNDU, "2.75406009586277478185e-01");
111 
112   check2("1.47116981099449883885", MPFR_RNDD, "3.86057874110010412760e-01");
113   check2("1.47116981099449883885", MPFR_RNDN, "3.86057874110010412760e-01");
114   check2("1.47116981099449883885", MPFR_RNDU, "3.86057874110010468272e-01");
115 
116   check2("1.58405446812987782401", MPFR_RNDD, "4.59987679246663727639e-01");
117   check2("1.58405446812987782401", MPFR_RNDN, "4.59987679246663783150e-01");
118   check2("1.58405446812987782401", MPFR_RNDU, "4.59987679246663783150e-01");
119 
120   check2("1.67192331263391547047", MPFR_RNDD, "5.13974647961076613889e-01");
121   check2("1.67192331263391547047", MPFR_RNDN, "5.13974647961076724911e-01");
122   check2("1.67192331263391547047", MPFR_RNDU, "5.13974647961076724911e-01");
123 
124   check2("1.71101198068990645318", MPFR_RNDD, "5.37084997042120315669e-01");
125   check2("1.71101198068990645318", MPFR_RNDN, "5.37084997042120315669e-01");
126   check2("1.71101198068990645318", MPFR_RNDU, "5.37084997042120426691e-01");
127 
128   check2("1.72634853551388700588", MPFR_RNDD, "5.46008504786553605648e-01");
129   check2("1.72634853551388700588", MPFR_RNDN, "5.46008504786553716670e-01");
130   check2("1.72634853551388700588", MPFR_RNDU, "5.46008504786553716670e-01");
131 
132   check2("2.00028876593004323325", MPFR_RNDD, "6.93291553102749702475e-01");
133   check2("2.00028876593004323325", MPFR_RNDN, "6.93291553102749813497e-01");
134   check2("2.00028876593004323325", MPFR_RNDU, "6.93291553102749813497e-01");
135 
136   check2("6.27593230200363105808", MPFR_RNDD, "1.83672204800630312072");
137   check2("6.27593230200363105808", MPFR_RNDN, "1.83672204800630334276");
138   check2("6.27593230200363105808", MPFR_RNDU, "1.83672204800630334276");
139 
140   check2("7.47216682321367997588", MPFR_RNDD, "2.01118502712453661729");
141   check2("7.47216682321367997588", MPFR_RNDN, "2.01118502712453706138");
142   check2("7.47216682321367997588", MPFR_RNDU, "2.01118502712453706138");
143 
144   check2("9.34589857718275318632", MPFR_RNDD, "2.23493759221664944903");
145   check2("9.34589857718275318632", MPFR_RNDN, "2.23493759221664989312");
146   check2("9.34589857718275318632", MPFR_RNDU, "2.23493759221664989312");
147 
148   check2("10.6856587560831854944", MPFR_RNDD, "2.36890253928838445674");
149   check2("10.6856587560831854944", MPFR_RNDN, "2.36890253928838445674");
150   check2("10.6856587560831854944", MPFR_RNDU, "2.36890253928838490083");
151 
152   check2("12.4646345033981766903", MPFR_RNDD, "2.52289539471636015122");
153   check2("12.4646345033981766903", MPFR_RNDN, "2.52289539471636015122");
154   check2("12.4646345033981766903", MPFR_RNDU, "2.52289539471636059531");
155 
156   check2("17.0953275851761752335", MPFR_RNDD, "2.83880518553861849185");
157   check2("17.0953275851761752335", MPFR_RNDN, "2.83880518553861893594");
158   check2("17.0953275851761752335", MPFR_RNDU, "2.83880518553861893594");
159 
160   check2("19.8509496207496916043", MPFR_RNDD, "2.98825184582516722998");
161   check2("19.8509496207496916043", MPFR_RNDN, "2.98825184582516722998");
162   check2("19.8509496207496916043", MPFR_RNDU, "2.98825184582516767406");
163 
164   check2("23.9512076062771335216", MPFR_RNDD, "3.17601874455977206679");
165   check2("23.9512076062771335216", MPFR_RNDN, "3.17601874455977206679");
166   check2("23.9512076062771335216", MPFR_RNDU, "3.17601874455977251088");
167 
168   check2("428.315247165198229595", MPFR_RNDD, "6.05985948325268264369");
169   check2("428.315247165198229595", MPFR_RNDN, "6.05985948325268353187");
170   check2("428.315247165198229595", MPFR_RNDU, "6.05985948325268353187");
171 }
172 
173 static void
special(void)174 special (void)
175 {
176   mpfr_t x, y;
177   int inex;
178   mpfr_exp_t emin, emax;
179   int r;
180 
181   emin = mpfr_get_emin ();
182   emax = mpfr_get_emax ();
183 
184   mpfr_init2 (x, 53);
185   mpfr_init2 (y, 53);
186 
187   /* Check special case: An overflow in const_pi could occurs! */
188   set_emin (-125);
189   set_emax (128);
190   mpfr_set_prec (y, 24*2);
191   mpfr_set_prec (x, 24);
192   mpfr_set_str_binary (x, "0.111110101010101011110101E0");
193   test_log (y, x, MPFR_RNDN);
194   set_emin (emin);
195   set_emax (emax);
196 
197   mpfr_set_prec (y, 53);
198   mpfr_set_prec (x, 53);
199   mpfr_set_ui (x, 3, MPFR_RNDD);
200   test_log (y, x, MPFR_RNDD);
201   if (mpfr_cmp_str1 (y, "1.09861228866810956"))
202     {
203       printf ("Error in mpfr_log(3) for MPFR_RNDD\n");
204       exit (1);
205     }
206 
207   /* check large precision */
208   mpfr_set_prec (x, 3322);
209   mpfr_set_prec (y, 3322);
210   mpfr_set_ui (x, 3, MPFR_RNDN);
211   mpfr_sqrt (x, x, MPFR_RNDN);
212   test_log (y, x, MPFR_RNDN);
213 
214   /* negative argument */
215   mpfr_set_si (x, -1, MPFR_RNDN);
216   test_log (y, x, MPFR_RNDN);
217   MPFR_ASSERTN(mpfr_nan_p (y));
218 
219   /* infinite loop when  */
220   set_emax (128);
221   mpfr_set_prec (x, 251);
222   mpfr_set_prec (y, 251);
223   mpfr_set_str_binary (x, "0.10010111000000000001101E8");
224   /* x = 4947981/32768, log(x) ~ 5.017282... */
225   test_log (y, x, MPFR_RNDN);
226 
227   set_emax (emax);
228 
229   mpfr_set_ui (x, 0, MPFR_RNDN);
230   inex = test_log (y, x, MPFR_RNDN);
231   MPFR_ASSERTN (inex == 0);
232   MPFR_ASSERTN (mpfr_inf_p (y));
233   MPFR_ASSERTN (mpfr_sgn (y) < 0);
234 
235   mpfr_set_ui (x, 0, MPFR_RNDN);
236   mpfr_neg (x, x, MPFR_RNDN);
237   inex = test_log (y, x, MPFR_RNDN);
238   MPFR_ASSERTN (inex == 0);
239   MPFR_ASSERTN (mpfr_inf_p (y));
240   MPFR_ASSERTN (mpfr_sgn (y) < 0);
241 
242   /* check log(1) is +0 whatever the rounding mode */
243   mpfr_set_ui (x, 1, MPFR_RNDN);
244   RND_LOOP (r)
245     {
246       mpfr_clear_flags ();
247       inex = test_log (y, x, (mpfr_rnd_t) r);
248       MPFR_ASSERTN (__gmpfr_flags == 0);
249       MPFR_ASSERTN (inex == 0);
250       MPFR_ASSERTN (MPFR_IS_ZERO (y));
251       MPFR_ASSERTN (MPFR_IS_POS (y));
252     }
253 
254   mpfr_clear (x);
255   mpfr_clear (y);
256 }
257 
258 static void
x_near_one(void)259 x_near_one (void)
260 {
261   mpfr_t x, y;
262   int inex;
263 
264   mpfr_init2 (x, 32);
265   mpfr_init2 (y, 16);
266 
267   mpfr_set_ui (x, 1, MPFR_RNDN);
268   mpfr_nextbelow (x);
269   inex = mpfr_log (y, x, MPFR_RNDD);
270   if (mpfr_cmp_str (y, "-0.1000000000000001E-31", 2, MPFR_RNDN)
271       || inex >= 0)
272     {
273       printf ("Failure in x_near_one, got inex = %d and\ny = ", inex);
274       mpfr_dump (y);
275     }
276 
277   mpfr_clears (x, y, (mpfr_ptr) 0);
278 }
279 
280 #define TEST_FUNCTION test_log
281 #define TEST_RANDOM_POS 8
282 #include "tgeneric.c"
283 
284 int
main(int argc,char * argv[])285 main (int argc, char *argv[])
286 {
287   tests_start_mpfr ();
288 
289   if (argc == 3 || argc == 4)
290     {   /* tlog x prec rnd */
291       check3 (argv[1], strtoul (argv[2], NULL, 10),
292               (argc == 4) ? (mpfr_rnd_t) atoi (argv[3]) : MPFR_RNDN);
293       goto done;
294     }
295 
296   special ();
297   check_worst_cases();
298 
299   check2("1.01979300812244555452", MPFR_RNDN, "1.95996734891603664741e-02");
300   check2("10.0",MPFR_RNDU,"2.30258509299404590110e+00");
301   check2("6.0",MPFR_RNDU,"1.79175946922805517936");
302   check2("1.0",MPFR_RNDZ,"0.0");
303   check2("62.0",MPFR_RNDU,"4.12713438504509166905");
304   check2("0.5",MPFR_RNDZ,"-6.93147180559945286226e-01");
305   check2("3.0",MPFR_RNDZ,"1.09861228866810956006e+00");
306   check2("234375765.0",MPFR_RNDU,"1.92724362186836231104e+01");
307   check2("8.0",MPFR_RNDZ,"2.07944154167983574765e+00");
308   check2("44.0",MPFR_RNDU,"3.78418963391826146392e+00");
309   check2("1.01979300812244555452", MPFR_RNDN, "1.95996734891603664741e-02");
310 
311   /* bugs found by Vincent Lefe`vre */
312   check2("0.99999599881598921769", MPFR_RNDN, "-0.0000040011920155404072924737977900999652547398000024259090423583984375");
313   check2("9.99995576063808955247e-01",MPFR_RNDZ,"-4.42394597667932383816e-06");
314   check2("9.99993687357856209097e-01",MPFR_RNDN,"-6.31266206860017342601e-06");
315   check2("9.99995223520736886691e-01",MPFR_RNDN,"-4.77649067052670982220e-06");
316   check2("9.99993025794720935551e-01",MPFR_RNDN,"-6.97422959894716163837e-06");
317   check2("9.99987549017837484833e-01",MPFR_RNDN,"-1.24510596766369924330e-05");
318   check2("9.99985901426543311032e-01",MPFR_RNDN,"-1.40986728425098585229e-05");
319   check2("9.99986053947420794330e-01",MPFR_RNDN, "-0.000013946149826301084938555592540598837558718514628708362579345703125");
320   check2("9.99971938247442126979e-01",MPFR_RNDN,"-2.80621462962173414790e-05");
321 
322   /* other bugs found by Vincent Lefe`vre */
323   check2("1.18615436389927785905e+77",MPFR_RNDN,"1.77469768607706015473e+02");
324   check2("9.48868723578399476187e+77",MPFR_RNDZ,"1.79549152432275803903e+02");
325   check2("2.31822210096938820854e+89",MPFR_RNDN,"2.05770873832573869322e+02");
326 
327   /* further bugs found by Vincent Lefe`vre */
328   check2("9.99999989485669482647e-01",MPFR_RNDZ,"-1.05143305726283042331e-08");
329   check2("9.99999989237970177136e-01",MPFR_RNDZ,"-1.07620298807745377934e-08");
330   check2("9.99999989239339082125e-01",MPFR_RNDN,"-1.07606609757704445430e-08");
331 
332   check2("7.3890560989306504",MPFR_RNDU,"2.0000000000000004"); /* exp(2.0) */
333   check2("7.3890560989306495",MPFR_RNDU,"2.0"); /* exp(2.0) */
334   check2("7.53428236571286402512e+34",MPFR_RNDZ,"8.03073567492226345621e+01");
335   check2("6.18784121531737948160e+19",MPFR_RNDZ,"4.55717030391710693493e+01");
336   check2("1.02560267603047283735e+00",MPFR_RNDD,"2.52804164149448735987e-02");
337   check2("7.53428236571286402512e+34",MPFR_RNDZ,"8.03073567492226345621e+01");
338   check2("1.42470900831881198052e+49",MPFR_RNDZ,"113.180637144887668910087086260318756103515625");
339 
340   check2("1.08013816255293777466e+11",MPFR_RNDN,"2.54055249841782604392e+01");
341   check2("6.72783635300509015581e-37",MPFR_RNDU,"-8.32893948416799503320e+01");
342   check2("2.25904918906057891180e-52",MPFR_RNDU,"-1.18919480823735682406e+02");
343   check2("1.48901209246462951085e+00",MPFR_RNDD,"3.98112874867437460668e-01");
344   check2("1.70322470467612341327e-01",MPFR_RNDN,"-1.77006175364294615626");
345   check2("1.94572026316065240791e+01",MPFR_RNDD,"2.96821731676437838842");
346   check2("4.01419512207026418764e+04",MPFR_RNDD,"1.06001772315501128218e+01");
347   check2("9.47077365236487591672e-04",MPFR_RNDZ,"-6.96212977303956748187e+00");
348   check2("3.95906157687589643802e-109",MPFR_RNDD,"-2.49605768114704119399e+02");
349   check2("2.73874914516503004113e-02",MPFR_RNDD,"-3.59766888618655977794e+00");
350   check2("9.18989072589566467669e-17",MPFR_RNDZ,"-3.69258425351464083519e+01");
351   check2("7706036453608191045959753324430048151991964994788917248.0",MPFR_RNDZ,"126.3815989984199177342816255986690521240234375");
352   check2("1.74827399630587801934e-23",MPFR_RNDZ,"-5.24008281254547156891e+01");
353   check2("4.35302958401482307665e+22",MPFR_RNDD,"5.21277441046519527390e+01");
354   check2("9.70791868689332915209e+00",MPFR_RNDD,"2.27294191194272210410e+00");
355   check2("2.22183639799464011100e-01",MPFR_RNDN,"-1.50425103275253957413e+00");
356   check2("2.27313466156682375540e+00",MPFR_RNDD,"8.21159787095675608448e-01");
357   check2("6.58057413965851156767e-01",MPFR_RNDZ,"-4.18463096196088235600e-01");
358   check2 ("7.34302197248998461006e+43",MPFR_RNDZ,"101.0049094695131799426235374994575977325439453125");
359   check2("6.09969788341579732815e+00",MPFR_RNDD,"1.80823924264386204363e+00");
360 
361   x_near_one ();
362 
363   test_generic (MPFR_PREC_MIN, 100, 40);
364 
365   data_check ("data/log", mpfr_log, "mpfr_log");
366   bad_cases (mpfr_log, mpfr_exp, "mpfr_log", 256, -30, 30, 4, 128, 800, 50);
367 
368  done:
369   tests_end_mpfr ();
370   return 0;
371 }
372