← Back to branch summary

~ubuntu-branches/debian/sid/genius/sid

« back to all changes in this revision

Viewing changes to mpfr/tests/tgamma.c

  • Committer: Bazaar Package Importer
  • Author(s): Daniel Holbach
  • Date: 2006-08-21 12:57:45 UTC
  • Revision ID: james.westby@ubuntu.com-20060821125745-sl9ks8v7fq324bdf
Tags: upstream-0.7.6.1
Import upstream version 0.7.6.1

Show diffs side-by-side

added added

removed removed

Lines of Context:
mpfr/tests/tgamma.c
 
1
/* mpfr_tgamma -- test file for gamma function
 
2
 
 
3
Copyright 2001, 2002, 2003, 2004, 2005 Free Software Foundation.
 
4
 
 
5
This file is part of the MPFR Library, and was contributed by Mathieu Dutour.
 
6
 
 
7
The MPFR Library is free software; you can redistribute it and/or modify
 
8
it under the terms of the GNU Lesser General Public License as published by
 
9
the Free Software Foundation; either version 2.1 of the License, or (at your
 
10
option) any later version.
 
11
 
 
12
The MPFR Library is distributed in the hope that it will be useful, but
 
13
WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
 
14
or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
 
15
License for more details.
 
16
 
 
17
You should have received a copy of the GNU Lesser General Public License
 
18
along with the MPFR Library; see the file COPYING.LIB.  If not, write to
 
19
the Free Software Foundation, Inc., 51 Franklin Place, Fifth Floor, Boston,
 
20
MA 02110-1301, USA. */
 
21
 
 
22
#include <stdio.h>
 
23
#include <stdlib.h>
 
24
 
 
25
#include "mpfr-test.h"
 
26
 
 
27
#define TEST_FUNCTION mpfr_gamma
 
28
#define RAND_FUNCTION(x) mpfr_random2(x, MPFR_LIMB_SIZE (x), 1)
 
29
#include "tgeneric.c"
 
30
 
 
31
static void
 
32
special (void)
 
33
{
 
34
  mpfr_t x, y;
 
35
  int inex;
 
36
 
 
37
  mpfr_init (x);
 
38
  mpfr_init (y);
 
39
 
 
40
  mpfr_set_nan (x);
 
41
  mpfr_gamma (y, x, GMP_RNDN);
 
42
  if (!mpfr_nan_p (y))
 
43
    {
 
44
      printf ("Error for gamma(NaN)\n");
 
45
      exit (1);
 
46
    }
 
47
 
 
48
  mpfr_set_inf (x, -1);
 
49
  mpfr_gamma (y, x, GMP_RNDN);
 
50
  if (!mpfr_nan_p (y))
 
51
    {
 
52
      printf ("Error for gamma(-Inf)\n");
 
53
      exit (1);
 
54
    }
 
55
 
 
56
  mpfr_set_inf (x, 1);
 
57
  mpfr_gamma (y, x, GMP_RNDN);
 
58
  if (!mpfr_inf_p (y) || mpfr_sgn (y) < 0)
 
59
    {
 
60
      printf ("Error for gamma(+Inf)\n");
 
61
      exit (1);
 
62
    }
 
63
 
 
64
  mpfr_set_ui (x, 0, GMP_RNDN);
 
65
  mpfr_gamma (y, x, GMP_RNDN);
 
66
  if (!mpfr_inf_p (y) || mpfr_sgn (y) < 0)
 
67
    {
 
68
      printf ("Error for gamma(+0)\n");
 
69
      exit (1);
 
70
    }
 
71
 
 
72
  mpfr_set_ui (x, 0, GMP_RNDN);
 
73
  mpfr_neg (x, x, GMP_RNDN);
 
74
  mpfr_gamma (y, x, GMP_RNDN);
 
75
  if (!mpfr_inf_p (y) || mpfr_sgn (y) > 0)
 
76
    {
 
77
      printf ("Error for gamma(-0)\n");
 
78
      exit (1);
 
79
    }
 
80
 
 
81
  mpfr_set_ui (x, 1, GMP_RNDN);
 
82
  mpfr_gamma (y, x, GMP_RNDN);
 
83
  if (mpfr_cmp_ui (y, 1))
 
84
    {
 
85
      printf ("Error for gamma(1)\n");
 
86
      exit (1);
 
87
    }
 
88
 
 
89
  mpfr_set_si (x, -1, GMP_RNDN);
 
90
  mpfr_gamma (y, x, GMP_RNDN);
 
91
  if (!mpfr_nan_p (y))
 
92
    {
 
93
      printf ("Error for gamma(-1)\n");
 
94
      exit (1);
 
95
    }
 
96
 
 
97
  mpfr_set_prec (x, 53);
 
98
  mpfr_set_prec (y, 53);
 
99
 
 
100
#define CHECK_X1 "1.0762904832837976166"
 
101
#define CHECK_Y1 "0.96134843256452096050"
 
102
 
 
103
  mpfr_set_str (x, CHECK_X1, 10, GMP_RNDN);
 
104
  mpfr_gamma (y, x, GMP_RNDN);
 
105
  mpfr_set_str (x, CHECK_Y1, 10, GMP_RNDN);
 
106
  if (mpfr_cmp (y, x))
 
107
    {
 
108
      printf ("mpfr_gamma("CHECK_X1") is wrong: expected %1.20e, got %1.20e\n",
 
109
              mpfr_get_d (x, GMP_RNDN), mpfr_get_d (y, GMP_RNDN));
 
110
      exit (1);
 
111
    }
 
112
 
 
113
#define CHECK_X2 "9.23709516716202383435e-01"
 
114
#define CHECK_Y2 "1.0502315560291053398"
 
115
  mpfr_set_str (x, CHECK_X2, 10, GMP_RNDN);
 
116
  mpfr_gamma (y, x, GMP_RNDN);
 
117
  mpfr_set_str (x, CHECK_Y2, 10, GMP_RNDN);
 
118
  if (mpfr_cmp (y, x))
 
119
    {
 
120
      printf ("mpfr_gamma("CHECK_X2") is wrong: expected %1.20e, got %1.20e\n",
 
121
              mpfr_get_d (x, GMP_RNDN), mpfr_get_d (y, GMP_RNDN));
 
122
      exit (1);
 
123
    }
 
124
 
 
125
  mpfr_set_prec (x, 8);
 
126
  mpfr_set_prec (y, 175);
 
127
  mpfr_set_ui (x, 33, GMP_RNDN);
 
128
  mpfr_gamma (y, x, GMP_RNDU);
 
129
  mpfr_set_prec (x, 175);
 
130
  mpfr_set_str_binary (x, "0.110010101011010101101000010101010111000110011101001000101011000001100010111001101001011E118");
 
131
  if (mpfr_cmp (x, y))
 
132
    {
 
133
      printf ("Error in mpfr_gamma (1)\n");
 
134
      exit (1);
 
135
    }
 
136
 
 
137
  mpfr_set_prec (x, 21);
 
138
  mpfr_set_prec (y, 8);
 
139
  mpfr_set_ui (y, 120, GMP_RNDN);
 
140
  mpfr_gamma (x, y, GMP_RNDZ);
 
141
  mpfr_set_prec (y, 21);
 
142
  mpfr_set_str_binary (y, "0.101111101110100110110E654");
 
143
  if (mpfr_cmp (x, y))
 
144
    {
 
145
      printf ("Error in mpfr_gamma (120)\n");
 
146
      printf ("Expected "); mpfr_print_binary (y); puts ("");
 
147
      printf ("Got      "); mpfr_print_binary (x); puts ("");
 
148
      exit (1);
 
149
    }
 
150
 
 
151
  mpfr_set_prec (x, 3);
 
152
  mpfr_set_prec (y, 206);
 
153
  mpfr_set_str_binary (x, "0.110e10");
 
154
  inex = mpfr_gamma (y, x, GMP_RNDN);
 
155
  mpfr_set_prec (x, 206);
 
156
  mpfr_set_str_binary (x, "0.110111100001000001101010010001000111000100000100111000010011100011011111001100011110101000111101101100110001001100110100001001111110000101010000100100011100010011101110000001000010001100010000101001111E6250");
 
157
  if (mpfr_cmp (x, y))
 
158
    {
 
159
      printf ("Error in mpfr_gamma (768)\n");
 
160
      exit (1);
 
161
    }
 
162
  if (inex <= 0)
 
163
    {
 
164
      printf ("Wrong flag for mpfr_gamma (768)\n");
 
165
      exit (1);
 
166
    }
 
167
 
 
168
  /* worst case to exercise retry */
 
169
  mpfr_set_prec (x, 1000);
 
170
  mpfr_set_prec (y, 869);
 
171
  mpfr_const_pi (x, GMP_RNDN);
 
172
  mpfr_gamma (y, x, GMP_RNDN);
 
173
 
 
174
  mpfr_set_prec (x, 4);
 
175
  mpfr_set_prec (y, 4);
 
176
  mpfr_set_str_binary (x, "-0.1100E-66");
 
177
  mpfr_gamma (y, x, GMP_RNDN);
 
178
  mpfr_set_str_binary (x, "-0.1011E67");
 
179
  if (mpfr_cmp (x, y))
 
180
    {
 
181
      printf ("Error for gamma(-0.1100E-66)\n");
 
182
      exit (1);
 
183
    }
 
184
 
 
185
  mpfr_clear (x);
 
186
  mpfr_clear (y);
 
187
}
 
188
 
 
189
static void
 
190
special_overflow (void)
 
191
{
 
192
  mpfr_t x, y;
 
193
  mp_exp_t emin = mpfr_get_emin ();
 
194
 
 
195
  set_emin (-125);
 
196
  set_emax (128);
 
197
 
 
198
  mpfr_init2 (x, 24);
 
199
  mpfr_init2 (y, 24);
 
200
  mpfr_set_str_binary (x, "0.101100100000000000110100E7");
 
201
  mpfr_gamma (y, x, GMP_RNDN);
 
202
  if (!mpfr_inf_p (y))
 
203
    {
 
204
      printf ("Overflow error.\n");
 
205
      mpfr_dump (y);
 
206
      exit (1);
 
207
    }
 
208
 
 
209
  /* problem mentioned by Kenneth Wilder, 18 Aug 2005 */
 
210
  mpfr_set_prec (x, 29);
 
211
  mpfr_set_prec (y, 29);
 
212
  mpfr_set_str (x, "-200000000.5", 10, GMP_RNDN); /* exact */
 
213
  mpfr_gamma (y, x, GMP_RNDN);
 
214
  if (!(mpfr_zero_p (y) && MPFR_SIGN (y) < 0))
 
215
    {
 
216
      printf ("Error for gamma(-200000000.5)\n");
 
217
      printf ("expected -0");
 
218
      printf ("got      ");
 
219
      mpfr_dump (y);
 
220
      exit (1);
 
221
    }
 
222
 
 
223
  mpfr_set_prec (x, 53);
 
224
  mpfr_set_prec (y, 53);
 
225
  mpfr_set_str (x, "-200000000.1", 10, GMP_RNDN);
 
226
  mpfr_gamma (y, x, GMP_RNDN);
 
227
  if (!(mpfr_zero_p (y) && MPFR_SIGN (y) < 0))
 
228
    {
 
229
      printf ("Error for gamma(-200000000.1), prec=53\n");
 
230
      printf ("expected -0");
 
231
      printf ("got      ");
 
232
      mpfr_dump (y);
 
233
      exit (1);
 
234
    }
 
235
 
 
236
  /* another problem mentioned by Kenneth Wilder, 29 Aug 2005 */
 
237
  mpfr_set_prec (x, 333);
 
238
  mpfr_set_prec (y, 14);
 
239
  mpfr_set_str (x, "-2.0000000000000000000000005", 10, GMP_RNDN);
 
240
  mpfr_gamma (y, x, GMP_RNDN);
 
241
  mpfr_set_prec (x, 14);
 
242
  mpfr_set_str_binary (x, "-11010011110001E66");
 
243
  if (mpfr_cmp (x, y))
 
244
    {
 
245
      printf ("Error for gamma(-2.0000000000000000000000005)\n");
 
246
      printf ("expected "); mpfr_dump (x);
 
247
      printf ("got      "); mpfr_dump (y);
 
248
      exit (1);
 
249
    }
 
250
 
 
251
  /* another tests from Kenneth Wilder, 31 Aug 2005 */
 
252
  set_emax (200);
 
253
  set_emin (-200);
 
254
  mpfr_set_prec (x, 38);
 
255
  mpfr_set_prec (y, 54);
 
256
  mpfr_set_str_binary (x, "0.11101111011100111101001001010110101001E-166");
 
257
  mpfr_gamma (y, x, GMP_RNDN);
 
258
  mpfr_set_prec (x, 54);
 
259
  mpfr_set_str_binary (x, "0.100010001101100001110110001010111111010000100101011E167");
 
260
  if (mpfr_cmp (x, y))
 
261
    {
 
262
      printf ("Error for gamma (test 1)\n");
 
263
      printf ("expected "); mpfr_dump (x);
 
264
      printf ("got      "); mpfr_dump (y);
 
265
      exit (1);
 
266
    }
 
267
 
 
268
  set_emax (1000);
 
269
  set_emin (-2000);
 
270
  mpfr_set_prec (x, 38);
 
271
  mpfr_set_prec (y, 71);
 
272
  mpfr_set_str_binary (x, "10101011011100001111111000010111110010E-1034");
 
273
  /* 184083777010*2^(-1034) */
 
274
  mpfr_gamma (y, x, GMP_RNDN);
 
275
  mpfr_set_prec (x, 71);
 
276
  mpfr_set_str_binary (x, "10111111001000011110010001000000000000110011110000000011101011111111100E926");
 
277
  /* 1762885132679550982140*2^926 */
 
278
  if (mpfr_cmp (x, y))
 
279
    {
 
280
      printf ("Error for gamma (test 2)\n");
 
281
      printf ("expected "); mpfr_dump (x);
 
282
      printf ("got      "); mpfr_dump (y);
 
283
      exit (1);
 
284
    }
 
285
 
 
286
  mpfr_set_prec (x, 38);
 
287
  mpfr_set_prec (y, 88);
 
288
  mpfr_set_str_binary (x, "10111100111001010000100001100100100101E-104");
 
289
  /* 202824096037*2^(-104) */
 
290
  mpfr_gamma (y, x, GMP_RNDN);
 
291
  mpfr_set_prec (x, 88);
 
292
  mpfr_set_str_binary (x, "1010110101111000111010111100010110101010100110111000001011000111000011101100001101110010E-21");
 
293
  /* 209715199999500283894743922*2^(-21) */
 
294
  if (mpfr_cmp (x, y))
 
295
    {
 
296
      printf ("Error for gamma (test 3)\n");
 
297
      printf ("expected "); mpfr_dump (x);
 
298
      printf ("got      "); mpfr_dump (y);
 
299
      exit (1);
 
300
    }
 
301
 
 
302
  mpfr_set_prec (x, 171);
 
303
  mpfr_set_prec (y, 38);
 
304
  mpfr_set_str (x, "-2993155353253689176481146537402947624254601559176535", 10,
 
305
                GMP_RNDN);
 
306
  mpfr_div_2exp (x, x, 170, GMP_RNDN);
 
307
  mpfr_gamma (y, x, GMP_RNDN);
 
308
  mpfr_set_prec (x, 38);
 
309
  mpfr_set_str (x, "201948391737", 10, GMP_RNDN);
 
310
  mpfr_mul_2exp (x, x, 92, GMP_RNDN);
 
311
  if (mpfr_cmp (x, y))
 
312
    {
 
313
      printf ("Error for gamma (test 5)\n");
 
314
      printf ("expected "); mpfr_dump (x);
 
315
      printf ("got      "); mpfr_dump (y);
 
316
      exit (1);
 
317
    }
 
318
 
 
319
  set_emin (-500000);
 
320
  mpfr_set_prec (x, 337);
 
321
  mpfr_set_prec (y, 38);
 
322
  mpfr_set_str (x, "-30000.000000000000000000000000000000000000000000001", 10,
 
323
                GMP_RNDN);
 
324
  mpfr_gamma (y, x, GMP_RNDN);
 
325
  mpfr_set_prec (x, 38);
 
326
  mpfr_set_str (x, "-3.623795987425E-121243", 10, GMP_RNDN);
 
327
  if (mpfr_cmp (x, y))
 
328
    {
 
329
      printf ("Error for gamma (test 7)\n");
 
330
      printf ("expected "); mpfr_dump (x);
 
331
      printf ("got      "); mpfr_dump (y);
 
332
      exit (1);
 
333
    }
 
334
 
 
335
  /* was producing infinite loop */
 
336
  set_emin (emin);
 
337
  mpfr_set_prec (x, 71);
 
338
  mpfr_set_prec (y, 71);
 
339
  mpfr_set_str (x, "-200000000.1", 10, GMP_RNDN);
 
340
  mpfr_gamma (y, x, GMP_RNDN);
 
341
  if (!(mpfr_zero_p (y) && MPFR_SIGN (y) < 0))
 
342
    {
 
343
      printf ("Error for gamma (test 8)\n");
 
344
      printf ("expected "); mpfr_dump (x);
 
345
      printf ("got      "); mpfr_dump (y);
 
346
      exit (1);
 
347
    }
 
348
 
 
349
  set_emax (1073741823);
 
350
  mpfr_set_prec (x, 29);
 
351
  mpfr_set_prec (y, 29);
 
352
  mpfr_set_str (x, "423786866", 10, GMP_RNDN);
 
353
  mpfr_gamma (y, x, GMP_RNDN);
 
354
  if (!mpfr_inf_p (y) || mpfr_sgn (y) < 0)
 
355
    {
 
356
      printf ("Error for gamma(423786866)\n");
 
357
      exit (1);
 
358
    }
 
359
 
 
360
  mpfr_clear (y);
 
361
  mpfr_clear (x);
 
362
  set_emin (MPFR_EMIN_MIN);
 
363
  set_emax (MPFR_EMAX_MAX);
 
364
}
 
365
 
 
366
int
 
367
main (void)
 
368
{
 
369
  MPFR_TEST_USE_RANDS ();
 
370
  tests_start_mpfr ();
 
371
 
 
372
  special ();
 
373
  special_overflow ();
 
374
  test_generic (2, 100, 2);
 
375
 
 
376
  tests_end_mpfr ();
 
377
  return 0;
 
378
}