[RTL-ifcvt] PR rtl-optimization/68506: Fix emitting order of insns in IF-THEN-JOIN...
[official-gcc.git] / libgfortran / io / write_float.def
blob183fe0b9f2ea683511bb29f9396be72141fe9905
1 /* Copyright (C) 2007-2015 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
6 This file is part of the GNU Fortran runtime library (libgfortran).
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.
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.
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.
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/>. */
27 #include "config.h"
29 typedef enum
30 { S_NONE, S_MINUS, S_PLUS }
31 sign_t;
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. */
36 static sign_t
37 calculate_sign (st_parameter_dt *dtp, int negative_flag)
39 sign_t s = S_NONE;
41 if (negative_flag)
42 s = S_MINUS;
43 else
44 switch (dtp->u.p.sign_status)
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;
58 return s;
62 /* Determine the precision except for EN format. For G format,
63 determines an upper bound to be used for sizing the buffer. */
65 static int
66 determine_precision (st_parameter_dt * dtp, const fnode * f, int len)
68 int precision = f->u.real.d;
70 switch (f->format)
72 case FMT_F:
73 case FMT_G:
74 precision += dtp->u.p.scale_factor;
75 break;
76 case FMT_ES:
77 /* Scale factor has no effect on output. */
78 break;
79 case FMT_E:
80 case FMT_D:
81 /* See F2008 10.7.2.3.3.6 */
82 if (dtp->u.p.scale_factor <= 0)
83 precision += dtp->u.p.scale_factor - 1;
84 break;
85 default:
86 return -1;
89 /* If the scale factor has a large negative value, we must do our
90 own rounding? Use ROUND='NEAREST', which should be what snprintf
91 is using as well. */
92 if (precision < 0 &&
93 (dtp->u.p.current_unit->round_status == ROUND_UNSPECIFIED
94 || dtp->u.p.current_unit->round_status == ROUND_PROCDEFINED))
95 dtp->u.p.current_unit->round_status = ROUND_NEAREST;
97 /* Add extra guard digits up to at least full precision when we do
98 our own rounding. */
99 if (dtp->u.p.current_unit->round_status != ROUND_UNSPECIFIED
100 && dtp->u.p.current_unit->round_status != ROUND_PROCDEFINED)
102 precision += 2 * len + 4;
103 if (precision < 0)
104 precision = 0;
107 return precision;
111 /* Output a real number according to its format which is FMT_G free. */
113 static bool
114 output_float (st_parameter_dt *dtp, const fnode *f, char *buffer, size_t size,
115 int nprinted, int precision, int sign_bit, bool zero_flag)
117 char *out;
118 char *digits;
119 int e, w, d, p, i;
120 char expchar, rchar;
121 format_token ft;
122 /* Number of digits before the decimal point. */
123 int nbefore;
124 /* Number of zeros after the decimal point. */
125 int nzero;
126 /* Number of digits after the decimal point. */
127 int nafter;
128 int leadzero;
129 int nblanks;
130 int ndigits, edigits;
131 sign_t sign;
133 ft = f->format;
134 w = f->u.real.w;
135 d = f->u.real.d;
136 p = dtp->u.p.scale_factor;
138 rchar = '5';
140 /* We should always know the field width and precision. */
141 if (d < 0)
142 internal_error (&dtp->common, "Unspecified precision");
144 sign = calculate_sign (dtp, sign_bit);
146 /* Calculate total number of digits. */
147 if (ft == FMT_F)
148 ndigits = nprinted - 2;
149 else
150 ndigits = precision + 1;
152 /* Read the exponent back in. */
153 if (ft != FMT_F)
154 e = atoi (&buffer[ndigits + 3]) + 1;
155 else
156 e = 0;
158 /* Make sure zero comes out as 0.0e0. */
159 if (zero_flag)
160 e = 0;
162 /* Normalize the fractional component. */
163 if (ft != FMT_F)
165 buffer[2] = buffer[1];
166 digits = &buffer[2];
168 else
169 digits = &buffer[1];
171 /* Figure out where to place the decimal point. */
172 switch (ft)
174 case FMT_F:
175 nbefore = ndigits - precision;
176 /* Make sure the decimal point is a '.'; depending on the
177 locale, this might not be the case otherwise. */
178 digits[nbefore] = '.';
179 if (p != 0)
181 if (p > 0)
184 memmove (digits + nbefore, digits + nbefore + 1, p);
185 digits[nbefore + p] = '.';
186 nbefore += p;
187 nafter = d - p;
188 if (nafter < 0)
189 nafter = 0;
190 nafter = d;
191 nzero = 0;
193 else /* p < 0 */
195 if (nbefore + p >= 0)
197 nzero = 0;
198 memmove (digits + nbefore + p + 1, digits + nbefore + p, -p);
199 nbefore += p;
200 digits[nbefore] = '.';
201 nafter = d;
203 else
205 nzero = -(nbefore + p);
206 memmove (digits + 1, digits, nbefore);
207 digits++;
208 nafter = d + nbefore;
209 nbefore = 0;
211 if (nzero > d)
212 nzero = d;
215 else
217 nzero = 0;
218 nafter = d;
221 while (digits[0] == '0' && nbefore > 0)
223 digits++;
224 nbefore--;
225 ndigits--;
228 expchar = 0;
229 /* If we need to do rounding ourselves, get rid of the dot by
230 moving the fractional part. */
231 if (dtp->u.p.current_unit->round_status != ROUND_UNSPECIFIED
232 && dtp->u.p.current_unit->round_status != ROUND_PROCDEFINED)
233 memmove (digits + nbefore, digits + nbefore + 1, ndigits - nbefore);
234 break;
236 case FMT_E:
237 case FMT_D:
238 i = dtp->u.p.scale_factor;
239 if (d <= 0 && p == 0)
241 generate_error (&dtp->common, LIBERROR_FORMAT, "Precision not "
242 "greater than zero in format specifier 'E' or 'D'");
243 return false;
245 if (p <= -d || p >= d + 2)
247 generate_error (&dtp->common, LIBERROR_FORMAT, "Scale factor "
248 "out of range in format specifier 'E' or 'D'");
249 return false;
252 if (!zero_flag)
253 e -= p;
254 if (p < 0)
256 nbefore = 0;
257 nzero = -p;
258 nafter = d + p;
260 else if (p > 0)
262 nbefore = p;
263 nzero = 0;
264 nafter = (d - p) + 1;
266 else /* p == 0 */
268 nbefore = 0;
269 nzero = 0;
270 nafter = d;
273 if (ft == FMT_E)
274 expchar = 'E';
275 else
276 expchar = 'D';
277 break;
279 case FMT_EN:
280 /* The exponent must be a multiple of three, with 1-3 digits before
281 the decimal point. */
282 if (!zero_flag)
283 e--;
284 if (e >= 0)
285 nbefore = e % 3;
286 else
288 nbefore = (-e) % 3;
289 if (nbefore != 0)
290 nbefore = 3 - nbefore;
292 e -= nbefore;
293 nbefore++;
294 nzero = 0;
295 nafter = d;
296 expchar = 'E';
297 break;
299 case FMT_ES:
300 if (!zero_flag)
301 e--;
302 nbefore = 1;
303 nzero = 0;
304 nafter = d;
305 expchar = 'E';
306 break;
308 default:
309 /* Should never happen. */
310 internal_error (&dtp->common, "Unexpected format token");
313 if (zero_flag)
314 goto skip;
316 /* Round the value. The value being rounded is an unsigned magnitude. */
317 switch (dtp->u.p.current_unit->round_status)
319 /* For processor defined and unspecified rounding we use
320 snprintf to print the exact number of digits needed, and thus
321 let snprintf handle the rounding. On system claiming support
322 for IEEE 754, this ought to be round to nearest, ties to
323 even, corresponding to the Fortran ROUND='NEAREST'. */
324 case ROUND_PROCDEFINED:
325 case ROUND_UNSPECIFIED:
326 case ROUND_ZERO: /* Do nothing and truncation occurs. */
327 goto skip;
328 case ROUND_UP:
329 if (sign_bit)
330 goto skip;
331 goto updown;
332 case ROUND_DOWN:
333 if (!sign_bit)
334 goto skip;
335 goto updown;
336 case ROUND_NEAREST:
337 /* Round compatible unless there is a tie. A tie is a 5 with
338 all trailing zero's. */
339 i = nafter + nbefore;
340 if (digits[i] == '5')
342 for(i++ ; i < ndigits; i++)
344 if (digits[i] != '0')
345 goto do_rnd;
347 /* It is a tie so round to even. */
348 switch (digits[nafter + nbefore - 1])
350 case '1':
351 case '3':
352 case '5':
353 case '7':
354 case '9':
355 /* If odd, round away from zero to even. */
356 break;
357 default:
358 /* If even, skip rounding, truncate to even. */
359 goto skip;
362 /* Fall through. */
363 /* The ROUND_COMPATIBLE is rounding away from zero when there is a tie. */
364 case ROUND_COMPATIBLE:
365 rchar = '5';
366 goto do_rnd;
369 updown:
371 rchar = '0';
372 if (ft != FMT_F && w > 0 && d == 0 && p == 0)
373 nbefore = 1;
374 /* Scan for trailing zeros to see if we really need to round it. */
375 for(i = nbefore + nafter; i < ndigits; i++)
377 if (digits[i] != '0')
378 goto do_rnd;
380 goto skip;
382 do_rnd:
384 if (nbefore + nafter == 0)
385 /* Handle the case Fw.0 and value < 1.0 */
387 ndigits = 0;
388 if (digits[0] >= rchar)
390 /* We rounded to zero but shouldn't have */
391 nbefore = 1;
392 digits--;
393 digits[0] = '1';
394 ndigits = 1;
397 else if (nbefore + nafter < ndigits)
399 i = ndigits = nbefore + nafter;
400 if (digits[i] >= rchar)
402 /* Propagate the carry. */
403 for (i--; i >= 0; i--)
405 if (digits[i] != '9')
407 digits[i]++;
408 break;
410 digits[i] = '0';
413 if (i < 0)
415 /* The carry overflowed. Fortunately we have some spare
416 space at the start of the buffer. We may discard some
417 digits, but this is ok because we already know they are
418 zero. */
419 digits--;
420 digits[0] = '1';
421 if (ft == FMT_F)
423 if (nzero > 0)
425 nzero--;
426 nafter++;
428 else
429 nbefore++;
431 else if (ft == FMT_EN)
433 nbefore++;
434 if (nbefore == 4)
436 nbefore = 1;
437 e += 3;
440 else
441 e++;
446 skip:
448 /* Calculate the format of the exponent field. */
449 if (expchar && !(dtp->u.p.g0_no_blanks && e == 0))
451 edigits = 1;
452 for (i = abs (e); i >= 10; i /= 10)
453 edigits++;
455 if (f->u.real.e < 0)
457 /* Width not specified. Must be no more than 3 digits. */
458 if (e > 999 || e < -999)
459 edigits = -1;
460 else
462 edigits = 4;
463 if (e > 99 || e < -99)
464 expchar = ' ';
467 else
469 /* Exponent width specified, check it is wide enough. */
470 if (edigits > f->u.real.e)
471 edigits = -1;
472 else
473 edigits = f->u.real.e + 2;
476 else
477 edigits = 0;
479 /* Scan the digits string and count the number of zeros. If we make it
480 all the way through the loop, we know the value is zero after the
481 rounding completed above. */
482 int hasdot = 0;
483 for (i = 0; i < ndigits + hasdot; i++)
485 if (digits[i] == '.')
486 hasdot = 1;
487 else if (digits[i] != '0')
488 break;
491 /* To format properly, we need to know if the rounded result is zero and if
492 so, we set the zero_flag which may have been already set for
493 actual zero. */
494 if (i == ndigits + hasdot)
496 zero_flag = true;
497 /* The output is zero, so set the sign according to the sign bit unless
498 -fno-sign-zero was specified. */
499 if (compile_options.sign_zero == 1)
500 sign = calculate_sign (dtp, sign_bit);
501 else
502 sign = calculate_sign (dtp, 0);
505 /* Pick a field size if none was specified, taking into account small
506 values that may have been rounded to zero. */
507 if (w <= 0)
509 if (zero_flag)
510 w = d + (sign != S_NONE ? 2 : 1) + (d == 0 ? 1 : 0);
511 else
513 w = nbefore + nzero + nafter + (sign != S_NONE ? 2 : 1);
514 w = w == 1 ? 2 : w;
518 /* Work out how much padding is needed. */
519 nblanks = w - (nbefore + nzero + nafter + edigits + 1);
520 if (sign != S_NONE)
521 nblanks--;
523 /* See if we have space for a zero before the decimal point. */
524 if (nbefore == 0 && nblanks > 0)
526 leadzero = 1;
527 nblanks--;
529 else
530 leadzero = 0;
532 if (dtp->u.p.g0_no_blanks)
534 w -= nblanks;
535 nblanks = 0;
538 /* Create the ouput buffer. */
539 out = write_block (dtp, w);
540 if (out == NULL)
541 return false;
543 /* Check the value fits in the specified field width. */
544 if (nblanks < 0 || edigits == -1 || w == 1 || (w == 2 && sign != S_NONE))
546 if (unlikely (is_char4_unit (dtp)))
548 gfc_char4_t *out4 = (gfc_char4_t *) out;
549 memset4 (out4, '*', w);
550 return false;
552 star_fill (out, w);
553 return false;
556 /* For internal character(kind=4) units, we duplicate the code used for
557 regular output slightly modified. This needs to be maintained
558 consistent with the regular code that follows this block. */
559 if (unlikely (is_char4_unit (dtp)))
561 gfc_char4_t *out4 = (gfc_char4_t *) out;
562 /* Pad to full field width. */
564 if ( ( nblanks > 0 ) && !dtp->u.p.no_leading_blank)
566 memset4 (out4, ' ', nblanks);
567 out4 += nblanks;
570 /* Output the initial sign (if any). */
571 if (sign == S_PLUS)
572 *(out4++) = '+';
573 else if (sign == S_MINUS)
574 *(out4++) = '-';
576 /* Output an optional leading zero. */
577 if (leadzero)
578 *(out4++) = '0';
580 /* Output the part before the decimal point, padding with zeros. */
581 if (nbefore > 0)
583 if (nbefore > ndigits)
585 i = ndigits;
586 memcpy4 (out4, digits, i);
587 ndigits = 0;
588 while (i < nbefore)
589 out4[i++] = '0';
591 else
593 i = nbefore;
594 memcpy4 (out4, digits, i);
595 ndigits -= i;
598 digits += i;
599 out4 += nbefore;
602 /* Output the decimal point. */
603 *(out4++) = dtp->u.p.current_unit->decimal_status
604 == DECIMAL_POINT ? '.' : ',';
605 if (ft == FMT_F
606 && (dtp->u.p.current_unit->round_status == ROUND_UNSPECIFIED
607 || dtp->u.p.current_unit->round_status == ROUND_PROCDEFINED))
608 digits++;
610 /* Output leading zeros after the decimal point. */
611 if (nzero > 0)
613 for (i = 0; i < nzero; i++)
614 *(out4++) = '0';
617 /* Output digits after the decimal point, padding with zeros. */
618 if (nafter > 0)
620 if (nafter > ndigits)
621 i = ndigits;
622 else
623 i = nafter;
625 memcpy4 (out4, digits, i);
626 while (i < nafter)
627 out4[i++] = '0';
629 digits += i;
630 ndigits -= i;
631 out4 += nafter;
634 /* Output the exponent. */
635 if (expchar && !(dtp->u.p.g0_no_blanks && e == 0))
637 if (expchar != ' ')
639 *(out4++) = expchar;
640 edigits--;
642 snprintf (buffer, size, "%+0*d", edigits, e);
643 memcpy4 (out4, buffer, edigits);
646 if (dtp->u.p.no_leading_blank)
648 out4 += edigits;
649 memset4 (out4, ' ' , nblanks);
650 dtp->u.p.no_leading_blank = 0;
652 return true;
653 } /* End of character(kind=4) internal unit code. */
655 /* Pad to full field width. */
657 if ( ( nblanks > 0 ) && !dtp->u.p.no_leading_blank)
659 memset (out, ' ', nblanks);
660 out += nblanks;
663 /* Output the initial sign (if any). */
664 if (sign == S_PLUS)
665 *(out++) = '+';
666 else if (sign == S_MINUS)
667 *(out++) = '-';
669 /* Output an optional leading zero. */
670 if (leadzero)
671 *(out++) = '0';
673 /* Output the part before the decimal point, padding with zeros. */
674 if (nbefore > 0)
676 if (nbefore > ndigits)
678 i = ndigits;
679 memcpy (out, digits, i);
680 ndigits = 0;
681 while (i < nbefore)
682 out[i++] = '0';
684 else
686 i = nbefore;
687 memcpy (out, digits, i);
688 ndigits -= i;
691 digits += i;
692 out += nbefore;
695 /* Output the decimal point. */
696 *(out++) = dtp->u.p.current_unit->decimal_status == DECIMAL_POINT ? '.' : ',';
697 if (ft == FMT_F
698 && (dtp->u.p.current_unit->round_status == ROUND_UNSPECIFIED
699 || dtp->u.p.current_unit->round_status == ROUND_PROCDEFINED))
700 digits++;
702 /* Output leading zeros after the decimal point. */
703 if (nzero > 0)
705 for (i = 0; i < nzero; i++)
706 *(out++) = '0';
709 /* Output digits after the decimal point, padding with zeros. */
710 if (nafter > 0)
712 if (nafter > ndigits)
713 i = ndigits;
714 else
715 i = nafter;
717 memcpy (out, digits, i);
718 while (i < nafter)
719 out[i++] = '0';
721 digits += i;
722 ndigits -= i;
723 out += nafter;
726 /* Output the exponent. */
727 if (expchar && !(dtp->u.p.g0_no_blanks && e == 0))
729 if (expchar != ' ')
731 *(out++) = expchar;
732 edigits--;
734 snprintf (buffer, size, "%+0*d", edigits, e);
735 memcpy (out, buffer, edigits);
738 if (dtp->u.p.no_leading_blank)
740 out += edigits;
741 memset( out , ' ' , nblanks );
742 dtp->u.p.no_leading_blank = 0;
745 return true;
749 /* Write "Infinite" or "Nan" as appropriate for the given format. */
751 static void
752 write_infnan (st_parameter_dt *dtp, const fnode *f, int isnan_flag, int sign_bit)
754 char * p, fin;
755 int nb = 0;
756 sign_t sign;
757 int mark;
759 if (f->format != FMT_B && f->format != FMT_O && f->format != FMT_Z)
761 sign = calculate_sign (dtp, sign_bit);
762 mark = (sign == S_PLUS || sign == S_MINUS) ? 8 : 7;
764 nb = f->u.real.w;
766 /* If the field width is zero, the processor must select a width
767 not zero. 4 is chosen to allow output of '-Inf' or '+Inf' */
769 if ((nb == 0) || dtp->u.p.g0_no_blanks)
771 if (isnan_flag)
772 nb = 3;
773 else
774 nb = (sign == S_PLUS || sign == S_MINUS) ? 4 : 3;
776 p = write_block (dtp, nb);
777 if (p == NULL)
778 return;
779 if (nb < 3)
781 if (unlikely (is_char4_unit (dtp)))
783 gfc_char4_t *p4 = (gfc_char4_t *) p;
784 memset4 (p4, '*', nb);
786 else
787 memset (p, '*', nb);
788 return;
791 if (unlikely (is_char4_unit (dtp)))
793 gfc_char4_t *p4 = (gfc_char4_t *) p;
794 memset4 (p4, ' ', nb);
796 else
797 memset(p, ' ', nb);
799 if (!isnan_flag)
801 if (sign_bit)
803 /* If the sign is negative and the width is 3, there is
804 insufficient room to output '-Inf', so output asterisks */
805 if (nb == 3)
807 if (unlikely (is_char4_unit (dtp)))
809 gfc_char4_t *p4 = (gfc_char4_t *) p;
810 memset4 (p4, '*', nb);
812 else
813 memset (p, '*', nb);
814 return;
816 /* The negative sign is mandatory */
817 fin = '-';
819 else
820 /* The positive sign is optional, but we output it for
821 consistency */
822 fin = '+';
824 if (unlikely (is_char4_unit (dtp)))
826 gfc_char4_t *p4 = (gfc_char4_t *) p;
828 if (nb > mark)
829 /* We have room, so output 'Infinity' */
830 memcpy4 (p4 + nb - 8, "Infinity", 8);
831 else
832 /* For the case of width equals mark, there is not enough room
833 for the sign and 'Infinity' so we go with 'Inf' */
834 memcpy4 (p4 + nb - 3, "Inf", 3);
836 if (sign == S_PLUS || sign == S_MINUS)
838 if (nb < 9 && nb > 3)
839 /* Put the sign in front of Inf */
840 p4[nb - 4] = (gfc_char4_t) fin;
841 else if (nb > 8)
842 /* Put the sign in front of Infinity */
843 p4[nb - 9] = (gfc_char4_t) fin;
845 return;
848 if (nb > mark)
849 /* We have room, so output 'Infinity' */
850 memcpy(p + nb - 8, "Infinity", 8);
851 else
852 /* For the case of width equals 8, there is not enough room
853 for the sign and 'Infinity' so we go with 'Inf' */
854 memcpy(p + nb - 3, "Inf", 3);
856 if (sign == S_PLUS || sign == S_MINUS)
858 if (nb < 9 && nb > 3)
859 p[nb - 4] = fin; /* Put the sign in front of Inf */
860 else if (nb > 8)
861 p[nb - 9] = fin; /* Put the sign in front of Infinity */
864 else
866 if (unlikely (is_char4_unit (dtp)))
868 gfc_char4_t *p4 = (gfc_char4_t *) p;
869 memcpy4 (p4 + nb - 3, "NaN", 3);
871 else
872 memcpy(p + nb - 3, "NaN", 3);
874 return;
879 /* Returns the value of 10**d. */
881 #define CALCULATE_EXP(x) \
882 static GFC_REAL_ ## x \
883 calculate_exp_ ## x (int d)\
885 int i;\
886 GFC_REAL_ ## x r = 1.0;\
887 for (i = 0; i< (d >= 0 ? d : -d); i++)\
888 r *= 10;\
889 r = (d >= 0) ? r : 1.0 / r;\
890 return r;\
893 CALCULATE_EXP(4)
895 CALCULATE_EXP(8)
897 #ifdef HAVE_GFC_REAL_10
898 CALCULATE_EXP(10)
899 #endif
901 #ifdef HAVE_GFC_REAL_16
902 CALCULATE_EXP(16)
903 #endif
904 #undef CALCULATE_EXP
907 /* Define a macro to build code for write_float. */
909 /* Note: Before output_float is called, snprintf is used to print to buffer the
910 number in the format +D.DDDDe+ddd.
912 # The result will always contain a decimal point, even if no
913 digits follow it
915 - The converted value is to be left adjusted on the field boundary
917 + A sign (+ or -) always be placed before a number
919 * prec is used as the precision
921 e format: [-]d.ddde±dd where there is one digit before the
922 decimal-point character and the number of digits after it is
923 equal to the precision. The exponent always contains at least two
924 digits; if the value is zero, the exponent is 00. */
927 #define TOKENPASTE(x, y) TOKENPASTE2(x, y)
928 #define TOKENPASTE2(x, y) x ## y
930 #define DTOA(suff,prec,val) TOKENPASTE(DTOA2,suff)(prec,val)
932 #define DTOA2(prec,val) \
933 snprintf (buffer, size, "%+-#.*e", (prec), (val))
935 #define DTOA2L(prec,val) \
936 snprintf (buffer, size, "%+-#.*Le", (prec), (val))
939 #if defined(GFC_REAL_16_IS_FLOAT128)
940 #define DTOA2Q(prec,val) \
941 __qmath_(quadmath_snprintf) (buffer, size, "%+-#.*Qe", (prec), (val))
942 #endif
944 #define FDTOA(suff,prec,val) TOKENPASTE(FDTOA2,suff)(prec,val)
946 /* For F format, we print to the buffer with f format. */
947 #define FDTOA2(prec,val) \
948 snprintf (buffer, size, "%+-#.*f", (prec), (val))
950 #define FDTOA2L(prec,val) \
951 snprintf (buffer, size, "%+-#.*Lf", (prec), (val))
954 #if defined(GFC_REAL_16_IS_FLOAT128)
955 #define FDTOA2Q(prec,val) \
956 __qmath_(quadmath_snprintf) (buffer, size, "%+-#.*Qf", \
957 (prec), (val))
958 #endif
962 /* Generate corresponding I/O format for FMT_G and output.
963 The rules to translate FMT_G to FMT_E or FMT_F from DEC fortran
964 LRM (table 11-2, Chapter 11, "I/O Formatting", P11-25) is:
966 Data Magnitude Equivalent Conversion
967 0< m < 0.1-0.5*10**(-d-1) Ew.d[Ee]
968 m = 0 F(w-n).(d-1), n' '
969 0.1-0.5*10**(-d-1)<= m < 1-0.5*10**(-d) F(w-n).d, n' '
970 1-0.5*10**(-d)<= m < 10-0.5*10**(-d+1) F(w-n).(d-1), n' '
971 10-0.5*10**(-d+1)<= m < 100-0.5*10**(-d+2) F(w-n).(d-2), n' '
972 ................ ..........
973 10**(d-1)-0.5*10**(-1)<= m <10**d-0.5 F(w-n).0,n(' ')
974 m >= 10**d-0.5 Ew.d[Ee]
976 notes: for Gw.d , n' ' means 4 blanks
977 for Gw.dEe, n' ' means e+2 blanks
978 for rounding modes adjustment, r, See Fortran F2008 10.7.5.2.2
979 the asm volatile is required for 32-bit x86 platforms. */
981 #define OUTPUT_FLOAT_FMT_G(x,y) \
982 static void \
983 output_float_FMT_G_ ## x (st_parameter_dt *dtp, const fnode *f, \
984 GFC_REAL_ ## x m, char *buffer, size_t size, \
985 int sign_bit, bool zero_flag, int comp_d) \
987 int e = f->u.real.e;\
988 int d = f->u.real.d;\
989 int w = f->u.real.w;\
990 fnode newf;\
991 GFC_REAL_ ## x exp_d, r = 0.5, r_sc;\
992 int low, high, mid;\
993 int ubound, lbound;\
994 char *p, pad = ' ';\
995 int save_scale_factor, nb = 0;\
996 bool result;\
997 int nprinted, precision;\
998 volatile GFC_REAL_ ## x temp;\
1000 save_scale_factor = dtp->u.p.scale_factor;\
1002 switch (dtp->u.p.current_unit->round_status)\
1004 case ROUND_ZERO:\
1005 r = sign_bit ? 1.0 : 0.0;\
1006 break;\
1007 case ROUND_UP:\
1008 r = 1.0;\
1009 break;\
1010 case ROUND_DOWN:\
1011 r = 0.0;\
1012 break;\
1013 default:\
1014 break;\
1017 exp_d = calculate_exp_ ## x (d);\
1018 r_sc = (1 - r / exp_d);\
1019 temp = 0.1 * r_sc;\
1020 if ((m > 0.0 && ((m < temp) || (r >= (exp_d - m))))\
1021 || ((m == 0.0) && !(compile_options.allow_std\
1022 & (GFC_STD_F2003 | GFC_STD_F2008)))\
1023 || d == 0)\
1025 newf.format = FMT_E;\
1026 newf.u.real.w = w;\
1027 newf.u.real.d = d - comp_d;\
1028 newf.u.real.e = e;\
1029 nb = 0;\
1030 precision = determine_precision (dtp, &newf, x);\
1031 nprinted = DTOA(y,precision,m); \
1032 goto finish;\
1035 mid = 0;\
1036 low = 0;\
1037 high = d + 1;\
1038 lbound = 0;\
1039 ubound = d + 1;\
1041 while (low <= high)\
1043 mid = (low + high) / 2;\
1045 temp = (calculate_exp_ ## x (mid - 1) * r_sc);\
1047 if (m < temp)\
1049 ubound = mid;\
1050 if (ubound == lbound + 1)\
1051 break;\
1052 high = mid - 1;\
1054 else if (m > temp)\
1056 lbound = mid;\
1057 if (ubound == lbound + 1)\
1059 mid ++;\
1060 break;\
1062 low = mid + 1;\
1064 else\
1066 mid++;\
1067 break;\
1071 nb = e <= 0 ? 4 : e + 2;\
1072 nb = nb >= w ? w - 1 : nb;\
1073 newf.format = FMT_F;\
1074 newf.u.real.w = w - nb;\
1075 newf.u.real.d = m == 0.0 ? d - 1 : -(mid - d - 1) ;\
1076 dtp->u.p.scale_factor = 0;\
1077 precision = determine_precision (dtp, &newf, x); \
1078 nprinted = FDTOA(y,precision,m); \
1080 finish:\
1081 result = output_float (dtp, &newf, buffer, size, nprinted, precision,\
1082 sign_bit, zero_flag);\
1083 dtp->u.p.scale_factor = save_scale_factor;\
1086 if (nb > 0 && !dtp->u.p.g0_no_blanks)\
1088 p = write_block (dtp, nb);\
1089 if (p == NULL)\
1090 return;\
1091 if (!result)\
1092 pad = '*';\
1093 if (unlikely (is_char4_unit (dtp)))\
1095 gfc_char4_t *p4 = (gfc_char4_t *) p;\
1096 memset4 (p4, pad, nb);\
1098 else \
1099 memset (p, pad, nb);\
1103 OUTPUT_FLOAT_FMT_G(4,)
1105 OUTPUT_FLOAT_FMT_G(8,)
1107 #ifdef HAVE_GFC_REAL_10
1108 OUTPUT_FLOAT_FMT_G(10,L)
1109 #endif
1111 #ifdef HAVE_GFC_REAL_16
1112 # ifdef GFC_REAL_16_IS_FLOAT128
1113 OUTPUT_FLOAT_FMT_G(16,Q)
1114 #else
1115 OUTPUT_FLOAT_FMT_G(16,L)
1116 #endif
1117 #endif
1119 #undef OUTPUT_FLOAT_FMT_G
1122 /* EN format is tricky since the number of significant digits depends
1123 on the magnitude. Solve it by first printing a temporary value and
1124 figure out the number of significant digits from the printed
1125 exponent. Values y, 0.95*10.0**e <= y <10.0**e, are rounded to
1126 10.0**e even when the final result will not be rounded to 10.0**e.
1127 For these values the exponent returned by atoi has to be decremented
1128 by one. The values y in the ranges
1129 (1000.0-0.5*10.0**(-d))*10.0**(3*n) <= y < 10.0*(3*(n+1))
1130 (100.0-0.5*10.0**(-d))*10.0**(3*n) <= y < 10.0*(3*n+2)
1131 (10.0-0.5*10.0**(-d))*10.0**(3*n) <= y < 10.0*(3*n+1)
1132 are correctly rounded respectively to 1.0...0*10.0*(3*(n+1)),
1133 100.0...0*10.0*(3*n), and 10.0...0*10.0*(3*n), where 0...0
1134 represents d zeroes, by the lines 279 to 297. */
1136 #define EN_PREC(x,y)\
1138 volatile GFC_REAL_ ## x tmp, one = 1.0;\
1139 tmp = * (GFC_REAL_ ## x *)source;\
1140 if (isfinite (tmp))\
1142 nprinted = DTOA(y,0,tmp);\
1143 int e = atoi (&buffer[4]);\
1144 if (buffer[1] == '1')\
1146 tmp = (calculate_exp_ ## x (-e)) * tmp;\
1147 tmp = one - (tmp < 0 ? -tmp : tmp); \
1148 if (tmp > 0)\
1149 e = e - 1;\
1151 nbefore = e%3;\
1152 if (nbefore < 0)\
1153 nbefore = 3 + nbefore;\
1155 else\
1156 nprinted = -1;\
1159 static int
1160 determine_en_precision (st_parameter_dt *dtp, const fnode *f,
1161 const char *source, int len)
1163 int nprinted;
1164 char buffer[10];
1165 const size_t size = 10;
1166 int nbefore; /* digits before decimal point - 1. */
1168 switch (len)
1170 case 4:
1171 EN_PREC(4,)
1172 break;
1174 case 8:
1175 EN_PREC(8,)
1176 break;
1178 #ifdef HAVE_GFC_REAL_10
1179 case 10:
1180 EN_PREC(10,L)
1181 break;
1182 #endif
1183 #ifdef HAVE_GFC_REAL_16
1184 case 16:
1185 # ifdef GFC_REAL_16_IS_FLOAT128
1186 EN_PREC(16,Q)
1187 # else
1188 EN_PREC(16,L)
1189 # endif
1190 break;
1191 #endif
1192 default:
1193 internal_error (NULL, "bad real kind");
1196 if (nprinted == -1)
1197 return -1;
1199 int prec = f->u.real.d + nbefore;
1200 if (dtp->u.p.current_unit->round_status != ROUND_UNSPECIFIED
1201 && dtp->u.p.current_unit->round_status != ROUND_PROCDEFINED)
1202 prec += 2 * len + 4;
1203 return prec;
1207 #define WRITE_FLOAT(x,y)\
1209 GFC_REAL_ ## x tmp;\
1210 tmp = * (GFC_REAL_ ## x *)source;\
1211 sign_bit = signbit (tmp);\
1212 if (!isfinite (tmp))\
1214 write_infnan (dtp, f, isnan (tmp), sign_bit);\
1215 return;\
1217 tmp = sign_bit ? -tmp : tmp;\
1218 zero_flag = (tmp == 0.0);\
1219 if (f->format == FMT_G)\
1220 output_float_FMT_G_ ## x (dtp, f, tmp, buffer, size, sign_bit, \
1221 zero_flag, comp_d);\
1222 else\
1224 if (f->format == FMT_F)\
1225 nprinted = FDTOA(y,precision,tmp); \
1226 else\
1227 nprinted = DTOA(y,precision,tmp); \
1228 output_float (dtp, f, buffer, size, nprinted, precision,\
1229 sign_bit, zero_flag);\
1233 /* Output a real number according to its format. */
1235 static void
1236 write_float (st_parameter_dt *dtp, const fnode *f, const char *source, \
1237 int len, int comp_d)
1239 int sign_bit, nprinted;
1240 int precision; /* Precision for snprintf call. */
1241 bool zero_flag;
1243 if (f->format != FMT_EN)
1244 precision = determine_precision (dtp, f, len);
1245 else
1246 precision = determine_en_precision (dtp, f, source, len);
1248 /* 4932 is the maximum exponent of long double and quad precision, 3
1249 extra characters for the sign, the decimal point, and the
1250 trailing null, and finally some extra digits depending on the
1251 requested precision. */
1252 const size_t size = 4932 + 3 + precision;
1253 #define BUF_STACK_SZ 5000
1254 char buf_stack[BUF_STACK_SZ];
1255 char *buffer;
1256 if (size > BUF_STACK_SZ)
1257 buffer = xmalloc (size);
1258 else
1259 buffer = buf_stack;
1261 switch (len)
1263 case 4:
1264 WRITE_FLOAT(4,)
1265 break;
1267 case 8:
1268 WRITE_FLOAT(8,)
1269 break;
1271 #ifdef HAVE_GFC_REAL_10
1272 case 10:
1273 WRITE_FLOAT(10,L)
1274 break;
1275 #endif
1276 #ifdef HAVE_GFC_REAL_16
1277 case 16:
1278 # ifdef GFC_REAL_16_IS_FLOAT128
1279 WRITE_FLOAT(16,Q)
1280 # else
1281 WRITE_FLOAT(16,L)
1282 # endif
1283 break;
1284 #endif
1285 default:
1286 internal_error (NULL, "bad real kind");
1288 if (size > BUF_STACK_SZ)
1289 free (buffer);