[BACK]Return to tpow3.c CVS log [TXT][DIR] Up to [local] / OpenXM_contrib / gmp / mpfr / tests

Annotation of OpenXM_contrib/gmp/mpfr/tests/tpow3.c, Revision 1.1.1.1

1.1       ohara       1: /* Test file for mpfr_pow.
                      2:
                      3: Copyright 2001 Free Software Foundation.
                      4: Adapted from tarctan.c.
                      5:
                      6: This file is part of the MPFR Library.
                      7:
                      8: The 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 2.1 of the License, or (at your
                     11: option) any later version.
                     12:
                     13: The 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 MPFR Library; see the file COPYING.LIB.  If not, write to
                     20: the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
                     21: MA 02111-1307, USA. */
                     22:
                     23: #include <stdio.h>
                     24: #include <limits.h>
                     25: #include <stdlib.h>
                     26: #include "gmp.h"
                     27: #include "gmp-impl.h"
                     28: #include "mpfr.h"
                     29: #include "mpfr-impl.h"
                     30: #include "mpfr-test.h"
                     31:
                     32:
                     33: int
                     34: main (int argc, char *argv[])
                     35: {
                     36:   mpfr_t x, y, z, ax;
                     37:   long int iy;
                     38:   mpfr_init (x);
                     39:   mpfr_init (ax);
                     40:   mpfr_init2 (y,sizeof(unsigned long int)*CHAR_BIT);
                     41:   mpfr_init (z);
                     42:
                     43:   MPFR_SET_NAN(x);
                     44:   mpfr_random(y);
                     45:   mpfr_pow (z, x,y, GMP_RNDN);
                     46:   if(!MPFR_IS_NAN(z))
                     47:     {
                     48:       printf ("evaluation of function in x=NAN does not return NAN");
                     49:       exit (1);
                     50:     }
                     51:
                     52:   MPFR_SET_NAN(y);
                     53:   mpfr_random(x);
                     54:   mpfr_pow (z, x,y, GMP_RNDN);
                     55:   if(!MPFR_IS_NAN(z))
                     56:     {
                     57:       printf ("evaluation of function in y=NAN does not return NAN");
                     58:       exit (1);
                     59:     }
                     60:
                     61:   MPFR_CLEAR_FLAGS(z);
                     62:   MPFR_CLEAR_FLAGS(y);
                     63:   MPFR_CLEAR_FLAGS(x);
                     64:
                     65:   MPFR_SET_ZERO(y);
                     66:   mpfr_random(x);
                     67:   mpfr_pow (z, x,y, GMP_RNDN);
                     68:   if(mpfr_cmp_ui(z,1)!=0 && !(MPFR_IS_NAN(x)))
                     69:     {
                     70:       printf ("evaluation of function in y=0 does not return 1\n");
                     71:       printf ("x =");
                     72:       mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
                     73:       printf ("\n y =");
                     74:       mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
                     75:       printf ("\n result =");
                     76:       mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
                     77:       exit (1);
                     78:     }
                     79:
                     80:   MPFR_CLEAR_FLAGS(z);
                     81:   MPFR_CLEAR_FLAGS(y);
                     82:   MPFR_CLEAR_FLAGS(x);
                     83:
                     84:   MPFR_SET_INF(y);
                     85:   if (MPFR_SIGN(y) < 0)
                     86:     MPFR_CHANGE_SIGN(y);
                     87:   mpfr_random(x);
                     88:   mpfr_set_prec (ax, MPFR_PREC(x));
                     89:   mpfr_abs(ax,x,GMP_RNDN);
                     90:   mpfr_pow (z, x,y, GMP_RNDN);
                     91:   if( !MPFR_IS_INF(z) && (mpfr_cmp_ui(ax,1) > 0) )
                     92:     {
                     93:       printf ("evaluation of function in y=INF (|x|>1) does not return INF");
                     94:       exit (1);
                     95:     }
                     96:   if( !MPFR_IS_ZERO(z) && (mpfr_cmp_ui(ax,1) < 0) )
                     97:     {
                     98:       printf ("\nevaluation of function in y=INF (|x|<1) does not return 0");
                     99:       printf ("\nx =");
                    100:       mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
                    101:       printf ("\n y =");
                    102:       mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
                    103:       printf ("\n result =");
                    104:       mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
                    105:       putchar('\n');
                    106:       exit (1);
                    107:     }
                    108:
                    109:
                    110:   MPFR_CLEAR_FLAGS(z);
                    111:   MPFR_CLEAR_FLAGS(y);
                    112:   MPFR_CLEAR_FLAGS(x);
                    113:
                    114:   MPFR_SET_INF(y);
                    115:   if (MPFR_SIGN(y) > 0)
                    116:     MPFR_CHANGE_SIGN(y);
                    117:   mpfr_random(x);
                    118:   mpfr_set_prec (ax, MPFR_PREC(x));
                    119:   mpfr_abs(ax,x,GMP_RNDN);
                    120:   mpfr_pow (z, x,y, GMP_RNDN);
                    121:   mpfr_pow (z, x,y, GMP_RNDN);
                    122:   if( !MPFR_IS_INF(z) && (mpfr_cmp_ui(ax,1) < 0) )
                    123:     {
                    124:       printf ("evaluation of function in y=INF (for |x| <0) does not return INF");
                    125:       exit (1);
                    126:     }
                    127:   if( !MPFR_IS_ZERO(z) && (mpfr_cmp_ui(ax,1) > 0) )
                    128:     {
                    129:       printf ("evaluation of function in y=INF (for |x| >0) does not return 0");
                    130:       exit (1);
                    131:     }
                    132:
                    133:   MPFR_CLEAR_FLAGS(z);
                    134:   MPFR_CLEAR_FLAGS(y);
                    135:   MPFR_CLEAR_FLAGS(x);
                    136:
                    137:   MPFR_SET_INF(x);
                    138:   if (MPFR_SIGN(x) < 0)
                    139:     MPFR_CHANGE_SIGN(x);
                    140:   mpfr_random(y);
                    141:   mpfr_pow (z, x,y, GMP_RNDN);
                    142:   if(!MPFR_IS_INF(z) && (MPFR_SIGN(y) > 0))
                    143:     {
                    144:       printf ("evaluation of function in INF does not return INF");
                    145:       printf ("\nx =");
                    146:       mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
                    147:       printf ("\n y =");
                    148:       mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
                    149:       printf ("\n result =");
                    150:       mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
                    151:       putchar('\n');
                    152:       exit (1);
                    153:     }
                    154:   if(!MPFR_IS_ZERO(z) && (MPFR_SIGN(y) < 0))
                    155:     {
                    156:       printf ("evaluation of function in INF does not return INF");
                    157:       printf ("\nx =");
                    158:       mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
                    159:       printf ("\n y =");
                    160:       mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
                    161:       printf ("\n result =");
                    162:       mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
                    163:       putchar('\n');
                    164:       exit (1);
                    165:     }
                    166:
                    167:
                    168:   MPFR_CLEAR_FLAGS(z);
                    169:   MPFR_CLEAR_FLAGS(y);
                    170:   MPFR_CLEAR_FLAGS(x);
                    171:
                    172:   MPFR_SET_INF(x);
                    173:   if (MPFR_SIGN(x) > 0)
                    174:     MPFR_CHANGE_SIGN(x);
                    175:   mpfr_random(y);
                    176:   if (random() % 2)
                    177:     mpfr_neg (y, y, GMP_RNDN);
                    178:    mpfr_pow (z, x,y, GMP_RNDN);
                    179:   if(!MPFR_IS_INF(z) && (MPFR_SIGN(y) > 0) && (mpfr_isinteger(y)))
                    180:     {
                    181:       printf ("evaluation of function in x=-INF does not return INF");
                    182:       printf ("\nx =");
                    183:       mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
                    184:       printf ("\n y =");
                    185:       mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
                    186:       printf ("\n result =");
                    187:       mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
                    188:       putchar('\n');
                    189:       if(mpfr_isinteger(y))
                    190:         printf("y is an integer\n");
                    191:       else
                    192:         printf("y is not an integer\n");
                    193:
                    194:       exit (1);
                    195:     }
                    196:   if(!MPFR_IS_ZERO(z) && (MPFR_SIGN(y) < 0) && (mpfr_isinteger(y)))
                    197:     {
                    198:       printf ("evaluation of function in x=-INF does not return 0");
                    199:       printf ("\nx =");
                    200:       mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
                    201:       printf ("\n y =");
                    202:       mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
                    203:       printf ("\n result =");
                    204:       mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
                    205:       putchar('\n');
                    206:
                    207:       if(mpfr_isinteger(y))
                    208:         printf("y is an integer\n");
                    209:       else
                    210:         printf("y is not an integer\n");
                    211:
                    212:       exit (1);
                    213:     }
                    214:   MPFR_CLEAR_FLAGS(z);
                    215:   MPFR_CLEAR_FLAGS(y);
                    216:   MPFR_CLEAR_FLAGS(x);
                    217:
                    218:   MPFR_SET_INF(x);
                    219:   if (MPFR_SIGN(x) > 0)
                    220:     MPFR_CHANGE_SIGN(x);
                    221:
                    222:   iy=random();
                    223:   mpfr_random(y);
                    224:   if (random() % 2)
                    225:     iy=-iy;
                    226:   mpfr_set_d(y,iy,GMP_RNDN);
                    227:   mpfr_pow (z, x,y, GMP_RNDN);
                    228:   if(!MPFR_IS_INF(z) && (MPFR_SIGN(y) > 0) && (mpfr_isinteger(y)))
                    229:     {
                    230:       printf ("evaluation of function in x=-INF does not return INF");
                    231:       printf ("\nx =");
                    232:       mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
                    233:       printf ("\n y =");
                    234:       mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
                    235:       printf ("\n result =");
                    236:       mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
                    237:       putchar('\n');
                    238:       if(mpfr_isinteger(y))
                    239:         printf("y is an integer\n");
                    240:       else
                    241:         printf("y is not an integer\n");
                    242:
                    243:       exit (1);
                    244:     }
                    245:   if(!MPFR_IS_ZERO(z) && (MPFR_SIGN(y) < 0) && (mpfr_isinteger(y)))
                    246:     {
                    247:       printf ("evaluation of function in x=-INF does not return 0");
                    248:       printf ("\nx =");
                    249:       mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
                    250:       printf ("\n y =");
                    251:       mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
                    252:       printf ("\n result =");
                    253:       mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
                    254:       putchar('\n');
                    255:
                    256:       if(mpfr_isinteger(y))
                    257:         printf("y is an integer\n");
                    258:       else
                    259:         printf("y is not an integer\n");
                    260:
                    261:       exit (1);
                    262:     }
                    263:   MPFR_CLEAR_FLAGS(z);
                    264:   MPFR_CLEAR_FLAGS(y);
                    265:   MPFR_CLEAR_FLAGS(x);
                    266:
                    267:   mpfr_set_ui(x,1,GMP_RNDN);
                    268:   MPFR_SET_INF(y);
                    269:   mpfr_pow (z, x,y, GMP_RNDN);
                    270:   if(!MPFR_IS_NAN(z))
                    271:     {
                    272:       printf ("evaluation of function in x=1, y=INF does not return NAN");
                    273:       printf ("\nx =");
                    274:       mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
                    275:       printf ("\n y =");
                    276:       mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
                    277:       printf ("\n result =");
                    278:       mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
                    279:       putchar('\n');
                    280:
                    281:       exit (1);
                    282:     }
                    283:
                    284:   MPFR_CLEAR_FLAGS(z);
                    285:   MPFR_CLEAR_FLAGS(y);
                    286:   MPFR_CLEAR_FLAGS(x);
                    287:
                    288:   MPFR_SET_ZERO(x);
                    289:   mpfr_random(y);
                    290:   if (random() % 2)
                    291:     mpfr_neg (y, y, GMP_RNDN);
                    292:
                    293:   mpfr_pow (z, x,y, GMP_RNDN);
                    294:
                    295:   if(!MPFR_IS_ZERO(z) && (MPFR_SIGN(y) < 0) && !(mpfr_isinteger(y)))
                    296:     {
                    297:       printf ("evaluation of function in y<0 does not return 0");
                    298:       printf ("\nx =");
                    299:       mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
                    300:       printf ("\n y =");
                    301:       mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
                    302:       printf ("\n result =");
                    303:       mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
                    304:       putchar('\n');
                    305:
                    306:       exit (1);
                    307:     }
                    308:   if(!MPFR_IS_INF(z) && (MPFR_SIGN(y) < 0) && (mpfr_isinteger(y)))
                    309:     {
                    310:       printf ("evaluation of function in y<0 (y integer) does not return INF");
                    311:       printf ("\nx =");
                    312:       mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
                    313:       printf ("\n y =");
                    314:       mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
                    315:       printf ("\n result =");
                    316:       mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
                    317:       putchar('\n');
                    318:       exit (1);
                    319:     }
                    320:   if(!MPFR_IS_ZERO(z) && (MPFR_SIGN(y) > 0) && (mpfr_isinteger(y)))
                    321:     {
                    322:       printf ("evaluation of function in y<0 (y integer) does not return 0");
                    323:        printf ("\nx =");
                    324:       mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
                    325:       printf ("\n y =");
                    326:       mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
                    327:       printf ("\n result =");
                    328:       mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
                    329:       putchar('\n');
                    330:      exit (1);
                    331:     }
                    332:
                    333:
                    334:   {
                    335:   mp_prec_t prec, yprec;
                    336:   mpfr_t t, s;
                    337:   mp_rnd_t rnd;
                    338:   int inexact, compare, compare2;
                    339:   unsigned int n, err;
                    340:
                    341:   int p0=2;
                    342:   int p1=100;
                    343:   int N=100;
                    344:
                    345:   mpfr_init (s);
                    346:   mpfr_init (t);
                    347:
                    348:   /* generic test */
                    349:   for (prec = p0; prec <= p1; prec++)
                    350:     {
                    351:       mpfr_set_prec (x, prec);
                    352:       mpfr_set_prec (s, sizeof(unsigned long int)*CHAR_BIT);
                    353:       mpfr_set_prec (z, prec);
                    354:       mpfr_set_prec (t, prec);
                    355:       yprec = prec + 10;
                    356:
                    357:       for (n=0; n<N; n++)
                    358:        {
                    359:
                    360:          mpfr_random (x);
                    361:
                    362:          mpfr_random (s);
                    363:           if (random() % 2)
                    364:             mpfr_neg (s, s, GMP_RNDN);
                    365:          rnd = random () % 4;
                    366:          mpfr_set_prec (y, yprec);
                    367:          compare = mpfr_pow (y, x, s, rnd);
                    368:          err = (rnd == GMP_RNDN) ? yprec + 1 : yprec;
                    369:          if (mpfr_can_round (y, err, rnd, rnd, prec))
                    370:            {
                    371:              mpfr_set (t, y, rnd);
                    372:              inexact = mpfr_pow (z,x, s, rnd);
                    373:              if (mpfr_cmp (t, z))
                    374:                {
                    375:                  printf ("results differ for x=");
                    376:                  mpfr_out_str (stdout, 2, prec, x, GMP_RNDN);
                    377:                   printf (" values of the exponential=");
                    378:                  mpfr_out_str (stdout, 2, prec, s, GMP_RNDN);
                    379:                  printf (" prec=%u rnd_mode=%s\n", (unsigned) prec,
                    380:                          mpfr_print_rnd_mode (rnd));
                    381:                  printf ("got      ");
                    382:                  mpfr_out_str (stdout, 2, prec, z, GMP_RNDN);
                    383:                  putchar ('\n');
                    384:                  printf ("expected ");
                    385:                  mpfr_out_str (stdout, 2, prec, t, GMP_RNDN);
                    386:                  putchar ('\n');
                    387:                  printf ("approx  ");
                    388:                  mpfr_print_binary (y);
                    389:                  putchar ('\n');
                    390:                  exit (1);
                    391:                }
                    392:              compare2 = mpfr_cmp (t, y);
                    393:              /* if rounding to nearest, cannot know the sign of t - f(x)
                    394:                 because of composed rounding: y = o(f(x)) and t = o(y) */
                    395:              if ((rnd != GMP_RNDN) && (compare * compare2 >= 0))
                    396:                compare = compare + compare2;
                    397:              else
                    398:                compare = inexact; /* cannot determine sign(t-f(x)) */
                    399:              if (((inexact == 0) && (compare != 0)) ||
                    400:                  ((inexact > 0) && (compare <= 0)) ||
                    401:                  ((inexact < 0) && (compare >= 0)))
                    402:                {
                    403:                  fprintf (stderr, "Wrong inexact flag for rnd=%s: expected %d, got %d\n",
                    404:                           mpfr_print_rnd_mode (rnd), compare, inexact);
                    405:                  printf ("x="); mpfr_print_binary (x); putchar ('\n');
                    406:                  printf ("y="); mpfr_print_binary (y); putchar ('\n');
                    407:                  printf ("t="); mpfr_print_binary (t); putchar ('\n');
                    408:                  exit (1);
                    409:                }
                    410:            }
                    411:        }
                    412:     }
                    413:
                    414:   mpfr_clear (s);
                    415:   mpfr_clear (t);
                    416:
                    417:   }
                    418:   mpfr_clear (x);
                    419:   mpfr_clear (y);
                    420:   mpfr_clear (z);
                    421:   mpfr_clear (ax);
                    422:
                    423:   return 0;
                    424: }

FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>