libgfortran/io/write_float.def

   1 /* Copyright (C) 2007, 2008, 2009, 2010, 2011 Free Software Foundation, Inc.
   2    Contributed by Andy Vaught
   3    Write float code factoring to this file by Jerry DeLisle
   4    F2003 I/O support contributed by Jerry DeLisle
   5
   6 This file is part of the GNU Fortran runtime library (libgfortran).
   7
   8 Libgfortran is free software; you can redistribute it and/or modify
   9 it under the terms of the GNU General Public License as published by
  10 the Free Software Foundation; either version 3, or (at your option)
  11 any later version.
  12
  13 Libgfortran is distributed in the hope that it will be useful,
  14 but WITHOUT ANY WARRANTY; without even the implied warranty of
  15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  16 GNU General Public License for more details.
  17
  18 Under Section 7 of GPL version 3, you are granted additional
  19 permissions described in the GCC Runtime Library Exception, version
  20 3.1, as published by the Free Software Foundation.
  21
  22 You should have received a copy of the GNU General Public License and
  23 a copy of the GCC Runtime Library Exception along with this program;
  24 see the files COPYING3 and COPYING.RUNTIME respectively.  If not, see
  25 <http://www.gnu.org/licenses/>.  */
  26
  27 #include "config.h"
  28
  29 typedef enum
  30 { S_NONE, S_MINUS, S_PLUS }
  31 sign_t;
  32
  33 /* Given a flag that indicates if a value is negative or not, return a
  34    sign_t that gives the sign that we need to produce.  */
  35
  36 static sign_t
  37 calculate_sign (st_parameter_dt *dtp, int negative_flag)
  38 {
  39   sign_t s = S_NONE;
  40
  41   if (negative_flag)
  42     s = S_MINUS;
  43   else
  44     switch (dtp->u.p.sign_status)
  45       {
  46       case SIGN_SP:     /* Show sign. */
  47         s = S_PLUS;
  48         break;
  49       case SIGN_SS:     /* Suppress sign. */
  50         s = S_NONE;
  51         break;
  52       case SIGN_S:      /* Processor defined. */
  53       case SIGN_UNSPECIFIED:
  54         s = options.optional_plus ? S_PLUS : S_NONE;
  55         break;
  56       }
  57
  58   return s;
  59 }
  60
  61
  62 /* Output a real number according to its format which is FMT_G free.  */
  63
  64 static try
  65 output_float (st_parameter_dt *dtp, const fnode *f, char *buffer, size_t size,
  66               int sign_bit, bool zero_flag, int ndigits, int edigits)
  67 {
  68   char *out;
  69   char *digits;
  70   int e;
  71   char expchar, rchar;
  72   format_token ft;
  73   int w;
  74   int d;
  75   /* Number of digits before the decimal point.  */
  76   int nbefore;
  77   /* Number of zeros after the decimal point.  */
  78   int nzero;
  79   /* Number of digits after the decimal point.  */
  80   int nafter;
  81   /* Number of zeros after the decimal point, whatever the precision.  */
  82   int nzero_real;
  83   int leadzero;
  84   int nblanks;
  85   int i;
  86   sign_t sign;
  87
  88   ft = f->format;
  89   w = f->u.real.w;
  90   d = f->u.real.d;
  91
  92   rchar = '5';
  93   nzero_real = -1;
  94
  95   /* We should always know the field width and precision.  */
  96   if (d < 0)
  97     internal_error (&dtp->common, "Unspecified precision");
  98
  99   sign = calculate_sign (dtp, sign_bit);
 100
 101   /* The following code checks the given string has punctuation in the correct
 102      places.  Uncomment if needed for debugging.
 103      if (d != 0 && ((buffer[2] != '.' && buffer[2] != ',')
 104                     || buffer[ndigits + 2] != 'e'))
 105        internal_error (&dtp->common, "printf is broken");  */
 106
 107   /* Read the exponent back in.  */
 108   e = atoi (&buffer[ndigits + 3]) + 1;
 109
 110   /* Make sure zero comes out as 0.0e0.   */
 111   if (zero_flag)
 112     {
 113       e = 0;
 114       if (compile_options.sign_zero != 1)
 115         sign = calculate_sign (dtp, 0);
 116
 117       /* Handle special cases.  */
 118       if (w == 0)
 119         w = d + (sign != S_NONE ? 2 : 1) + (d == 0 ? 1 : 0);
 120     }
 121
 122   /* Normalize the fractional component.  */
 123   buffer[2] = buffer[1];
 124   digits = &buffer[2];
 125
 126   /* Figure out where to place the decimal point.  */
 127   switch (ft)
 128     {
 129     case FMT_F:
 130       if (d == 0 && e <= 0 && dtp->u.p.scale_factor == 0)
 131         {
 132           memmove (digits + 1, digits, ndigits - 1);
 133           digits[0] = '0';
 134           e++;
 135         }
 136
 137       nbefore = e + dtp->u.p.scale_factor;
 138       if (nbefore < 0)
 139         {
 140           nzero = -nbefore;
 141           nzero_real = nzero;
 142           if (nzero > d)
 143             nzero = d;
 144           nafter = d - nzero;
 145           nbefore = 0;
 146         }
 147       else
 148         {
 149           nzero = 0;
 150           nafter = d;
 151         }
 152       expchar = 0;
 153       break;
 154
 155     case FMT_E:
 156     case FMT_D:
 157       i = dtp->u.p.scale_factor;
 158       if (d <= 0 && i == 0)
 159         {
 160           generate_error (&dtp->common, LIBERROR_FORMAT, "Precision not "
 161                           "greater than zero in format specifier 'E' or 'D'");
 162           return FAILURE;
 163         }
 164       if (i <= -d || i >= d + 2)
 165         {
 166           generate_error (&dtp->common, LIBERROR_FORMAT, "Scale factor "
 167                           "out of range in format specifier 'E' or 'D'");
 168           return FAILURE;
 169         }
 170
 171       if (!zero_flag)
 172         e -= i;
 173       if (i < 0)
 174         {
 175           nbefore = 0;
 176           nzero = -i;
 177           nafter = d + i;
 178         }
 179       else if (i > 0)
 180         {
 181           nbefore = i;
 182           nzero = 0;
 183           nafter = (d - i) + 1;
 184         }
 185       else /* i == 0 */
 186         {
 187           nbefore = 0;
 188           nzero = 0;
 189           nafter = d;
 190         }
 191
 192       if (ft == FMT_E)
 193         expchar = 'E';
 194       else
 195         expchar = 'D';
 196       break;
 197
 198     case FMT_EN:
 199       /* The exponent must be a multiple of three, with 1-3 digits before
 200          the decimal point.  */
 201       if (!zero_flag)
 202         e--;
 203       if (e >= 0)
 204         nbefore = e % 3;
 205       else
 206         {
 207           nbefore = (-e) % 3;
 208           if (nbefore != 0)
 209             nbefore = 3 - nbefore;
 210         }
 211       e -= nbefore;
 212       nbefore++;
 213       nzero = 0;
 214       nafter = d;
 215       expchar = 'E';
 216       break;
 217
 218     case FMT_ES:
 219       if (!zero_flag)
 220         e--;
 221       nbefore = 1;
 222       nzero = 0;
 223       nafter = d;
 224       expchar = 'E';
 225       break;
 226
 227     default:
 228       /* Should never happen.  */
 229       internal_error (&dtp->common, "Unexpected format token");
 230     }
 231
 232   /* Round the value.  The value being rounded is an unsigned magnitude.
 233      The ROUND_COMPATIBLE is rounding away from zero when there is a tie.  */
 234   switch (dtp->u.p.current_unit->round_status)
 235     {
 236       case ROUND_ZERO: /* Do nothing and truncation occurs.  */
 237         goto skip;
 238       case ROUND_UP:
 239         if (sign_bit)
 240           goto skip;
 241         rchar = '0';
 242         break;
 243       case ROUND_DOWN:
 244         if (!sign_bit)
 245           goto skip;
 246         rchar = '0';
 247         break;
 248       case ROUND_NEAREST:
 249         /* Round compatible unless there is a tie. A tie is a 5 with
 250            all trailing zero's.  */
 251         i = nafter + nbefore;
 252         if (digits[i] == '5')
 253           {
 254             for(i++ ; i < ndigits; i++)
 255               {
 256                 if (digits[i] != '0')
 257                   goto do_rnd;
 258               }
 259             /* It is a  tie so round to even.  */
 260             switch (digits[nafter + nbefore - 1])
 261               {
 262                 case '1':
 263                 case '3':
 264                 case '5':
 265                 case '7':
 266                 case '9':
 267                   /* If odd, round away from zero to even.  */
 268                   break;
 269                 default:
 270                   /* If even, skip rounding, truncate to even.  */
 271                   goto skip;
 272               }
 273           }
 274          /* Fall through.  */
 275       case ROUND_PROCDEFINED:
 276       case ROUND_UNSPECIFIED:
 277       case ROUND_COMPATIBLE:
 278         rchar = '5';
 279         /* Just fall through and do the actual rounding.  */
 280     }
 281
 282   do_rnd:
 283
 284   if (nbefore + nafter == 0)
 285     {
 286       ndigits = 0;
 287       if (nzero_real == d && digits[0] >= rchar)
 288         {
 289           /* We rounded to zero but shouldn't have */
 290           nzero--;
 291           nafter = 1;
 292           digits[0] = '1';
 293           ndigits = 1;
 294         }
 295     }
 296   else if (nbefore + nafter < ndigits)
 297     {
 298       ndigits = nbefore + nafter;
 299       i = ndigits;
 300       if (digits[i] >= rchar)
 301         {
 302           /* Propagate the carry.  */
 303           for (i--; i >= 0; i--)
 304             {
 305               if (digits[i] != '9')
 306                 {
 307                   digits[i]++;
 308                   break;
 309                 }
 310               digits[i] = '0';
 311             }
 312
 313           if (i < 0)
 314             {
 315               /* The carry overflowed.  Fortunately we have some spare
 316                  space at the start of the buffer.  We may discard some
 317                  digits, but this is ok because we already know they are
 318                  zero.  */
 319               digits--;
 320               digits[0] = '1';
 321               if (ft == FMT_F)
 322                 {
 323                   if (nzero > 0)
 324                     {
 325                       nzero--;
 326                       nafter++;
 327                     }
 328                   else
 329                     nbefore++;
 330                 }
 331               else if (ft == FMT_EN)
 332                 {
 333                   nbefore++;
 334                   if (nbefore == 4)
 335                     {
 336                       nbefore = 1;
 337                       e += 3;
 338                     }
 339                 }
 340               else
 341                 e++;
 342             }
 343         }
 344     }
 345
 346   skip:
 347
 348   /* Calculate the format of the exponent field.  */
 349   if (expchar)
 350     {
 351       edigits = 1;
 352       for (i = abs (e); i >= 10; i /= 10)
 353         edigits++;
 354
 355       if (f->u.real.e < 0)
 356         {
 357           /* Width not specified.  Must be no more than 3 digits.  */
 358           if (e > 999 || e < -999)
 359             edigits = -1;
 360           else
 361             {
 362               edigits = 4;
 363               if (e > 99 || e < -99)
 364                 expchar = ' ';
 365             }
 366         }
 367       else
 368         {
 369           /* Exponent width specified, check it is wide enough.  */
 370           if (edigits > f->u.real.e)
 371             edigits = -1;
 372           else
 373             edigits = f->u.real.e + 2;
 374         }
 375     }
 376   else
 377     edigits = 0;
 378
 379   /* Zero values always output as positive, even if the value was negative
 380      before rounding.  */
 381   for (i = 0; i < ndigits; i++)
 382     {
 383       if (digits[i] != '0')
 384         break;
 385     }
 386   if (i == ndigits)
 387     {
 388       /* The output is zero, so set the sign according to the sign bit unless
 389          -fno-sign-zero was specified.  */
 390       if (compile_options.sign_zero == 1)
 391         sign = calculate_sign (dtp, sign_bit);
 392       else
 393         sign = calculate_sign (dtp, 0);
 394     }
 395
 396   /* Pick a field size if none was specified.  */
 397   if (w <= 0)
 398     {
 399       w = nbefore + nzero + nafter + (sign != S_NONE ? 2 : 1);
 400       w = w == 1 ? 2 : w;
 401     }
 402
 403   /* Work out how much padding is needed.  */
 404   nblanks = w - (nbefore + nzero + nafter + edigits + 1);
 405   if (sign != S_NONE)
 406     nblanks--;
 407
 408   if (dtp->u.p.g0_no_blanks)
 409     {
 410       w -= nblanks;
 411       nblanks = 0;
 412     }
 413
 414   /* Create the ouput buffer.  */
 415   out = write_block (dtp, w);
 416   if (out == NULL)
 417     return FAILURE;
 418
 419   /* Check the value fits in the specified field width.  */
 420   if (nblanks < 0 || edigits == -1 || w == 1 || (w == 2 && sign != S_NONE))
 421     {
 422       if (unlikely (is_char4_unit (dtp)))
 423         {
 424           gfc_char4_t *out4 = (gfc_char4_t *) out;
 425           memset4 (out4, '*', w);
 426           return FAILURE;
 427         }
 428       star_fill (out, w);
 429       return FAILURE;
 430     }
 431
 432   /* See if we have space for a zero before the decimal point.  */
 433   if (nbefore == 0 && nblanks > 0)
 434     {
 435       leadzero = 1;
 436       nblanks--;
 437     }
 438   else
 439     leadzero = 0;
 440
 441   /* For internal character(kind=4) units, we duplicate the code used for
 442      regular output slightly modified.  This needs to be maintained
 443      consistent with the regular code that follows this block.  */
 444   if (unlikely (is_char4_unit (dtp)))
 445     {
 446       gfc_char4_t *out4 = (gfc_char4_t *) out;
 447       /* Pad to full field width.  */
 448
 449       if ( ( nblanks > 0 ) && !dtp->u.p.no_leading_blank)
 450         {
 451           memset4 (out4, ' ', nblanks);
 452           out4 += nblanks;
 453         }
 454
 455       /* Output the initial sign (if any).  */
 456       if (sign == S_PLUS)
 457         *(out4++) = '+';
 458       else if (sign == S_MINUS)
 459         *(out4++) = '-';
 460
 461       /* Output an optional leading zero.  */
 462       if (leadzero)
 463         *(out4++) = '0';
 464
 465       /* Output the part before the decimal point, padding with zeros.  */
 466       if (nbefore > 0)
 467         {
 468           if (nbefore > ndigits)
 469             {
 470               i = ndigits;
 471               memcpy4 (out4, digits, i);
 472               ndigits = 0;
 473               while (i < nbefore)
 474                 out4[i++] = '0';
 475             }
 476           else
 477             {
 478               i = nbefore;
 479               memcpy4 (out4, digits, i);
 480               ndigits -= i;
 481             }
 482
 483           digits += i;
 484           out4 += nbefore;
 485         }
 486
 487       /* Output the decimal point.  */
 488       *(out4++) = dtp->u.p.current_unit->decimal_status
 489                     == DECIMAL_POINT ? '.' : ',';
 490
 491       /* Output leading zeros after the decimal point.  */
 492       if (nzero > 0)
 493         {
 494           for (i = 0; i < nzero; i++)
 495             *(out4++) = '0';
 496         }
 497
 498       /* Output digits after the decimal point, padding with zeros.  */
 499       if (nafter > 0)
 500         {
 501           if (nafter > ndigits)
 502             i = ndigits;
 503           else
 504             i = nafter;
 505
 506           memcpy4 (out4, digits, i);
 507           while (i < nafter)
 508             out4[i++] = '0';
 509
 510           digits += i;
 511           ndigits -= i;
 512           out4 += nafter;
 513         }
 514
 515       /* Output the exponent.  */
 516       if (expchar)
 517         {
 518           if (expchar != ' ')
 519             {
 520               *(out4++) = expchar;
 521               edigits--;
 522             }
 523 #if HAVE_SNPRINTF
 524           snprintf (buffer, size, "%+0*d", edigits, e);
 525 #else
 526           sprintf (buffer, "%+0*d", edigits, e);
 527 #endif
 528           memcpy4 (out4, buffer, edigits);
 529         }
 530
 531       if (dtp->u.p.no_leading_blank)
 532         {
 533           out4 += edigits;
 534           memset4 (out4, ' ' , nblanks);
 535           dtp->u.p.no_leading_blank = 0;
 536         }
 537       return SUCCESS;
 538     } /* End of character(kind=4) internal unit code.  */
 539
 540   /* Pad to full field width.  */
 541
 542   if ( ( nblanks > 0 ) && !dtp->u.p.no_leading_blank)
 543     {
 544       memset (out, ' ', nblanks);
 545       out += nblanks;
 546     }
 547
 548   /* Output the initial sign (if any).  */
 549   if (sign == S_PLUS)
 550     *(out++) = '+';
 551   else if (sign == S_MINUS)
 552     *(out++) = '-';
 553
 554   /* Output an optional leading zero.  */
 555   if (leadzero)
 556     *(out++) = '0';
 557
 558   /* Output the part before the decimal point, padding with zeros.  */
 559   if (nbefore > 0)
 560     {
 561       if (nbefore > ndigits)
 562         {
 563           i = ndigits;
 564           memcpy (out, digits, i);
 565           ndigits = 0;
 566           while (i < nbefore)
 567             out[i++] = '0';
 568         }
 569       else
 570         {
 571           i = nbefore;
 572           memcpy (out, digits, i);
 573           ndigits -= i;
 574         }
 575
 576       digits += i;
 577       out += nbefore;
 578     }
 579
 580   /* Output the decimal point.  */
 581   *(out++) = dtp->u.p.current_unit->decimal_status == DECIMAL_POINT ? '.' : ',';
 582
 583   /* Output leading zeros after the decimal point.  */
 584   if (nzero > 0)
 585     {
 586       for (i = 0; i < nzero; i++)
 587         *(out++) = '0';
 588     }
 589
 590   /* Output digits after the decimal point, padding with zeros.  */
 591   if (nafter > 0)
 592     {
 593       if (nafter > ndigits)
 594         i = ndigits;
 595       else
 596         i = nafter;
 597
 598       memcpy (out, digits, i);
 599       while (i < nafter)
 600         out[i++] = '0';
 601
 602       digits += i;
 603       ndigits -= i;
 604       out += nafter;
 605     }
 606
 607   /* Output the exponent.  */
 608   if (expchar)
 609     {
 610       if (expchar != ' ')
 611         {
 612           *(out++) = expchar;
 613           edigits--;
 614         }
 615 #if HAVE_SNPRINTF
 616       snprintf (buffer, size, "%+0*d", edigits, e);
 617 #else
 618       sprintf (buffer, "%+0*d", edigits, e);
 619 #endif
 620       memcpy (out, buffer, edigits);
 621     }
 622
 623   if (dtp->u.p.no_leading_blank)
 624     {
 625       out += edigits;
 626       memset( out , ' ' , nblanks );
 627       dtp->u.p.no_leading_blank = 0;
 628     }
 629
 630 #undef STR
 631 #undef STR1
 632 #undef MIN_FIELD_WIDTH
 633   return SUCCESS;
 634 }
 635
 636
 637 /* Write "Infinite" or "Nan" as appropriate for the given format.  */
 638
 639 static void
 640 write_infnan (st_parameter_dt *dtp, const fnode *f, int isnan_flag, int sign_bit)
 641 {
 642   char * p, fin;
 643   int nb = 0;
 644   sign_t sign;
 645   int mark;
 646
 647   if (f->format != FMT_B && f->format != FMT_O && f->format != FMT_Z)
 648     {
 649       sign = calculate_sign (dtp, sign_bit);
 650       mark = (sign == S_PLUS || sign == S_MINUS) ? 8 : 7;
 651
 652       nb =  f->u.real.w;
 653
 654       /* If the field width is zero, the processor must select a width
 655          not zero.  4 is chosen to allow output of '-Inf' or '+Inf' */
 656
 657       if (nb == 0)
 658         {
 659           if (isnan_flag)
 660             nb = 3;
 661           else
 662             nb = (sign == S_PLUS || sign == S_MINUS) ? 4 : 3;
 663         }
 664       p = write_block (dtp, nb);
 665       if (p == NULL)
 666         return;
 667       if (nb < 3)
 668         {
 669           if (unlikely (is_char4_unit (dtp)))
 670             {
 671               gfc_char4_t *p4 = (gfc_char4_t *) p;
 672               memset4 (p4, '*', nb);
 673             }
 674           else
 675             memset (p, '*', nb);
 676           return;
 677         }
 678
 679       if (unlikely (is_char4_unit (dtp)))
 680         {
 681           gfc_char4_t *p4 = (gfc_char4_t *) p;
 682           memset4 (p4, ' ', nb);
 683         }
 684       else
 685         memset(p, ' ', nb);
 686
 687       if (!isnan_flag)
 688         {
 689           if (sign_bit)
 690             {
 691               /* If the sign is negative and the width is 3, there is
 692                  insufficient room to output '-Inf', so output asterisks */
 693               if (nb == 3)
 694                 {
 695                   if (unlikely (is_char4_unit (dtp)))
 696                     {
 697                       gfc_char4_t *p4 = (gfc_char4_t *) p;
 698                       memset4 (p4, '*', nb);
 699                     }
 700                   else
 701                     memset (p, '*', nb);
 702                   return;
 703                 }
 704               /* The negative sign is mandatory */
 705               fin = '-';
 706             }
 707           else
 708             /* The positive sign is optional, but we output it for
 709                consistency */
 710             fin = '+';
 711
 712           if (unlikely (is_char4_unit (dtp)))
 713             {
 714               gfc_char4_t *p4 = (gfc_char4_t *) p;
 715
 716               if (nb > mark)
 717                 /* We have room, so output 'Infinity' */
 718                 memcpy4 (p4 + nb - 8, "Infinity", 8);
 719               else
 720                 /* For the case of width equals mark, there is not enough room
 721                    for the sign and 'Infinity' so we go with 'Inf' */
 722                 memcpy4 (p4 + nb - 3, "Inf", 3);
 723
 724               if (sign == S_PLUS || sign == S_MINUS)
 725                 {
 726                   if (nb < 9 && nb > 3)
 727                     /* Put the sign in front of Inf */
 728                     p4[nb - 4] = (gfc_char4_t) fin;
 729                   else if (nb > 8)
 730                     /* Put the sign in front of Infinity */
 731                     p4[nb - 9] = (gfc_char4_t) fin;
 732                 }
 733               return;
 734             }
 735
 736           if (nb > mark)
 737             /* We have room, so output 'Infinity' */
 738             memcpy(p + nb - 8, "Infinity", 8);
 739           else
 740             /* For the case of width equals 8, there is not enough room
 741                for the sign and 'Infinity' so we go with 'Inf' */
 742             memcpy(p + nb - 3, "Inf", 3);
 743
 744           if (sign == S_PLUS || sign == S_MINUS)
 745             {
 746               if (nb < 9 && nb > 3)
 747                 p[nb - 4] = fin;  /* Put the sign in front of Inf */
 748               else if (nb > 8)
 749                 p[nb - 9] = fin;  /* Put the sign in front of Infinity */
 750             }
 751         }
 752       else
 753         {
 754           if (unlikely (is_char4_unit (dtp)))
 755             {
 756               gfc_char4_t *p4 = (gfc_char4_t *) p;
 757               memcpy4 (p4 + nb - 3, "NaN", 3);
 758             }
 759           else
 760             memcpy(p + nb - 3, "NaN", 3);
 761         }
 762       return;
 763     }
 764 }
 765
 766
 767 /* Returns the value of 10**d.  */
 768
 769 #define CALCULATE_EXP(x) \
 770 inline static GFC_REAL_ ## x \
 771 calculate_exp_ ## x  (int d)\
 772 {\
 773   int i;\
 774   GFC_REAL_ ## x r = 1.0;\
 775   for (i = 0; i< (d >= 0 ? d : -d); i++)\
 776     r *= 10;\
 777   r = (d >= 0) ? r : 1.0 / r;\
 778   return r;\
 779 }
 780
 781 CALCULATE_EXP(4)
 782
 783 CALCULATE_EXP(8)
 784
 785 #ifdef HAVE_GFC_REAL_10
 786 CALCULATE_EXP(10)
 787 #endif
 788
 789 #ifdef HAVE_GFC_REAL_16
 790 CALCULATE_EXP(16)
 791 #endif
 792 #undef CALCULATE_EXP
 793
 794 /* Generate corresponding I/O format for FMT_G and output.
 795    The rules to translate FMT_G to FMT_E or FMT_F from DEC fortran
 796    LRM (table 11-2, Chapter 11, "I/O Formatting", P11-25) is:
 797
 798    Data Magnitude                              Equivalent Conversion
 799    0< m < 0.1-0.5*10**(-d-1)                   Ew.d[Ee]
 800    m = 0                                       F(w-n).(d-1), n' '
 801    0.1-0.5*10**(-d-1)<= m < 1-0.5*10**(-d)     F(w-n).d, n' '
 802    1-0.5*10**(-d)<= m < 10-0.5*10**(-d+1)      F(w-n).(d-1), n' '
 803    10-0.5*10**(-d+1)<= m < 100-0.5*10**(-d+2)  F(w-n).(d-2), n' '
 804    ................                           ..........
 805    10**(d-1)-0.5*10**(-1)<= m <10**d-0.5       F(w-n).0,n(' ')
 806    m >= 10**d-0.5                              Ew.d[Ee]
 807
 808    notes: for Gw.d ,  n' ' means 4 blanks
 809           for Gw.dEe, n' ' means e+2 blanks  */
 810
 811 #define OUTPUT_FLOAT_FMT_G(x) \
 812 static void \
 813 output_float_FMT_G_ ## x (st_parameter_dt *dtp, const fnode *f, \
 814                       GFC_REAL_ ## x m, char *buffer, size_t size, \
 815                       int sign_bit, bool zero_flag, int ndigits, int edigits) \
 816 { \
 817   int e = f->u.real.e;\
 818   int d = f->u.real.d;\
 819   int w = f->u.real.w;\
 820   fnode *newf;\
 821   GFC_REAL_ ## x rexp_d;\
 822   int low, high, mid;\
 823   int ubound, lbound;\
 824   char *p, pad = ' ';\
 825   int save_scale_factor, nb = 0;\
 826   try result;\
 827 \
 828   save_scale_factor = dtp->u.p.scale_factor;\
 829   newf = (fnode *) get_mem (sizeof (fnode));\
 830 \
 831   rexp_d = calculate_exp_ ## x (-d);\
 832   if ((m > 0.0 && m < 0.1 - 0.05 * rexp_d) || (rexp_d * (m + 0.5) >= 1.0) ||\
 833       ((m == 0.0) && !(compile_options.allow_std & GFC_STD_F2003)))\
 834     { \
 835       newf->format = FMT_E;\
 836       newf->u.real.w = w;\
 837       newf->u.real.d = d;\
 838       newf->u.real.e = e;\
 839       nb = 0;\
 840       goto finish;\
 841     }\
 842 \
 843   mid = 0;\
 844   low = 0;\
 845   high = d + 1;\
 846   lbound = 0;\
 847   ubound = d + 1;\
 848 \
 849   while (low <= high)\
 850     { \
 851       GFC_REAL_ ## x temp;\
 852       mid = (low + high) / 2;\
 853 \
 854       temp = (calculate_exp_ ## x (mid - 1) * (1 - 0.5 * rexp_d));\
 855 \
 856       if (m < temp)\
 857         { \
 858           ubound = mid;\
 859           if (ubound == lbound + 1)\
 860             break;\
 861           high = mid - 1;\
 862         }\
 863       else if (m > temp)\
 864         { \
 865           lbound = mid;\
 866           if (ubound == lbound + 1)\
 867             { \
 868               mid ++;\
 869               break;\
 870             }\
 871           low = mid + 1;\
 872         }\
 873       else\
 874         {\
 875           mid++;\
 876           break;\
 877         }\
 878     }\
 879 \
 880   if (e > 4)\
 881     e = 4;\
 882   if (e < 0)\
 883     nb = 4;\
 884   else\
 885     nb = e + 2;\
 886 \
 887   nb = nb >= w ? 0 : nb;\
 888   newf->format = FMT_F;\
 889   newf->u.real.w = f->u.real.w - nb;\
 890 \
 891   if (m == 0.0)\
 892     newf->u.real.d = d - 1;\
 893   else\
 894     newf->u.real.d = - (mid - d - 1);\
 895 \
 896   dtp->u.p.scale_factor = 0;\
 897 \
 898  finish:\
 899   result = output_float (dtp, newf, buffer, size, sign_bit, zero_flag, \
 900                          ndigits, edigits);\
 901   dtp->u.p.scale_factor = save_scale_factor;\
 902 \
 903   free (newf);\
 904 \
 905   if (nb > 0 && !dtp->u.p.g0_no_blanks)\
 906     {\
 907       p = write_block (dtp, nb);\
 908       if (p == NULL)\
 909         return;\
 910       if (result == FAILURE)\
 911         pad = '*';\
 912       if (unlikely (is_char4_unit (dtp)))\
 913         {\
 914           gfc_char4_t *p4 = (gfc_char4_t *) p;\
 915           memset4 (p4, pad, nb);\
 916         }\
 917       else\
 918         memset (p, pad, nb);\
 919     }\
 920 }\
 921
 922 OUTPUT_FLOAT_FMT_G(4)
 923
 924 OUTPUT_FLOAT_FMT_G(8)
 925
 926 #ifdef HAVE_GFC_REAL_10
 927 OUTPUT_FLOAT_FMT_G(10)
 928 #endif
 929
 930 #ifdef HAVE_GFC_REAL_16
 931 OUTPUT_FLOAT_FMT_G(16)
 932 #endif
 933
 934 #undef OUTPUT_FLOAT_FMT_G
 935
 936
 937 /* Define a macro to build code for write_float.  */
 938
 939   /* Note: Before output_float is called, sprintf is used to print to buffer the
 940      number in the format +D.DDDDe+ddd. For an N digit exponent, this gives us
 941      (MIN_FIELD_WIDTH-5)-N digits after the decimal point, plus another one
 942      before the decimal point.
 943
 944      #   The result will always contain a decimal point, even if no
 945          digits follow it
 946
 947      -   The converted value is to be left adjusted on the field boundary
 948
 949      +   A sign (+ or -) always be placed before a number
 950
 951      MIN_FIELD_WIDTH  minimum field width
 952
 953      *   (ndigits-1) is used as the precision
 954
 955      e format: [-]d.dddeÂ±dd where there is one digit before the
 956        decimal-point character and the number of digits after it is
 957        equal to the precision. The exponent always contains at least two
 958        digits; if the value is zero, the exponent is 00.  */
 959
 960 #ifdef HAVE_SNPRINTF
 961
 962 #define DTOA \
 963 snprintf (buffer, size, "%+-#" STR(MIN_FIELD_WIDTH) ".*" \
 964           "e", ndigits - 1, tmp);
 965
 966 #define DTOAL \
 967 snprintf (buffer, size, "%+-#" STR(MIN_FIELD_WIDTH) ".*" \
 968           "Le", ndigits - 1, tmp);
 969
 970 #else
 971
 972 #define DTOA \
 973 sprintf (buffer, "%+-#" STR(MIN_FIELD_WIDTH) ".*" \
 974          "e", ndigits - 1, tmp);
 975
 976 #define DTOAL \
 977 sprintf (buffer, "%+-#" STR(MIN_FIELD_WIDTH) ".*" \
 978          "Le", ndigits - 1, tmp);
 979
 980 #endif
 981
 982 #if defined(GFC_REAL_16_IS_FLOAT128)
 983 #define DTOAQ \
 984 __qmath_(quadmath_snprintf) (buffer, sizeof buffer, \
 985                              "%+-#" STR(MIN_FIELD_WIDTH) ".*" \
 986                              "Qe", ndigits - 1, tmp);
 987 #endif
 988
 989 #define WRITE_FLOAT(x,y)\
 990 {\
 991         GFC_REAL_ ## x tmp;\
 992         tmp = * (GFC_REAL_ ## x *)source;\
 993         sign_bit = signbit (tmp);\
 994         if (!isfinite (tmp))\
 995           { \
 996             write_infnan (dtp, f, isnan (tmp), sign_bit);\
 997             return;\
 998           }\
 999         tmp = sign_bit ? -tmp : tmp;\
1000         zero_flag = (tmp == 0.0);\
1001 \
1002         DTOA ## y\
1003 \
1004         if (f->format != FMT_G)\
1005           output_float (dtp, f, buffer, size, sign_bit, zero_flag, ndigits, \
1006                         edigits);\
1007         else \
1008           output_float_FMT_G_ ## x (dtp, f, tmp, buffer, size, sign_bit, \
1009                                     zero_flag, ndigits, edigits);\
1010 }\
1011
1012 /* Output a real number according to its format.  */
1013
1014 static void
1015 write_float (st_parameter_dt *dtp, const fnode *f, const char *source, int len)
1016 {
1017
1018 #if defined(HAVE_GFC_REAL_16) || __LDBL_DIG__ > 18
1019 # define MIN_FIELD_WIDTH 46
1020 #else
1021 # define MIN_FIELD_WIDTH 31
1022 #endif
1023 #define STR(x) STR1(x)
1024 #define STR1(x) #x
1025
1026   /* This must be large enough to accurately hold any value.  */
1027   char buffer[MIN_FIELD_WIDTH+1];
1028   int sign_bit, ndigits, edigits;
1029   bool zero_flag;
1030   size_t size;
1031
1032   size = MIN_FIELD_WIDTH+1;
1033
1034   /* printf pads blanks for us on the exponent so we just need it big enough
1035      to handle the largest number of exponent digits expected.  */
1036   edigits=4;
1037
1038   if (f->format == FMT_F || f->format == FMT_EN || f->format == FMT_G
1039       || ((f->format == FMT_D || f->format == FMT_E)
1040       && dtp->u.p.scale_factor != 0))
1041     {
1042       /* Always convert at full precision to avoid double rounding.  */
1043       ndigits = MIN_FIELD_WIDTH - 4 - edigits;
1044     }
1045   else
1046     {
1047       /* The number of digits is known, so let printf do the rounding.  */
1048       if (f->format == FMT_ES)
1049         ndigits = f->u.real.d + 1;
1050       else
1051         ndigits = f->u.real.d;
1052       if (ndigits > MIN_FIELD_WIDTH - 4 - edigits)
1053         ndigits = MIN_FIELD_WIDTH - 4 - edigits;
1054     }
1055
1056   switch (len)
1057     {
1058     case 4:
1059       WRITE_FLOAT(4,)
1060       break;
1061
1062     case 8:
1063       WRITE_FLOAT(8,)
1064       break;
1065
1066 #ifdef HAVE_GFC_REAL_10
1067     case 10:
1068       WRITE_FLOAT(10,L)
1069       break;
1070 #endif
1071 #ifdef HAVE_GFC_REAL_16
1072     case 16:
1073 # ifdef GFC_REAL_16_IS_FLOAT128
1074       WRITE_FLOAT(16,Q)
1075 # else
1076       WRITE_FLOAT(16,L)
1077 # endif
1078       break;
1079 #endif
1080     default:
1081       internal_error (NULL, "bad real kind");
1082     }
1083 }