-/* Copyright (C) 2007, 2008, 2009, 2010 Free Software Foundation, Inc.
+/* Copyright (C) 2007-2018 Free Software Foundation, Inc.
Contributed by Andy Vaught
Write float code factoring to this file by Jerry DeLisle
F2003 I/O support contributed by Jerry DeLisle
}
-/* Output a real number according to its format which is FMT_G free. */
+/* Determine the precision except for EN format. For G format,
+ determines an upper bound to be used for sizing the buffer. */
+
+static int
+determine_precision (st_parameter_dt * dtp, const fnode * f, int len)
+{
+ int precision = f->u.real.d;
+
+ switch (f->format)
+ {
+ case FMT_F:
+ case FMT_G:
+ precision += dtp->u.p.scale_factor;
+ break;
+ case FMT_ES:
+ /* Scale factor has no effect on output. */
+ break;
+ case FMT_E:
+ case FMT_D:
+ /* See F2008 10.7.2.3.3.6 */
+ if (dtp->u.p.scale_factor <= 0)
+ precision += dtp->u.p.scale_factor - 1;
+ break;
+ default:
+ return -1;
+ }
+
+ /* If the scale factor has a large negative value, we must do our
+ own rounding? Use ROUND='NEAREST', which should be what snprintf
+ is using as well. */
+ if (precision < 0 &&
+ (dtp->u.p.current_unit->round_status == ROUND_UNSPECIFIED
+ || dtp->u.p.current_unit->round_status == ROUND_PROCDEFINED))
+ dtp->u.p.current_unit->round_status = ROUND_NEAREST;
+
+ /* Add extra guard digits up to at least full precision when we do
+ our own rounding. */
+ if (dtp->u.p.current_unit->round_status != ROUND_UNSPECIFIED
+ && dtp->u.p.current_unit->round_status != ROUND_PROCDEFINED)
+ {
+ precision += 2 * len + 4;
+ if (precision < 0)
+ precision = 0;
+ }
+
+ return precision;
+}
+
+
+/* Build a real number according to its format which is FMT_G free. */
static void
-output_float (st_parameter_dt *dtp, const fnode *f, char *buffer, size_t size,
- int sign_bit, bool zero_flag, int ndigits, int edigits)
+build_float_string (st_parameter_dt *dtp, const fnode *f, char *buffer,
+ size_t size, int nprinted, int precision, int sign_bit,
+ bool zero_flag, int npad, char *result, size_t *len)
{
- char *out;
+ char *put;
char *digits;
- int e;
+ int e, w, d, p, i;
char expchar, rchar;
format_token ft;
- int w;
- int d;
/* Number of digits before the decimal point. */
int nbefore;
/* Number of zeros after the decimal point. */
int nzero;
/* Number of digits after the decimal point. */
int nafter;
- /* Number of zeros after the decimal point, whatever the precision. */
- int nzero_real;
int leadzero;
int nblanks;
- int i;
+ int ndigits, edigits;
sign_t sign;
ft = f->format;
w = f->u.real.w;
d = f->u.real.d;
+ p = dtp->u.p.scale_factor;
rchar = '5';
- nzero_real = -1;
/* We should always know the field width and precision. */
if (d < 0)
sign = calculate_sign (dtp, sign_bit);
- /* The following code checks the given string has punctuation in the correct
- places. Uncomment if needed for debugging.
- if (d != 0 && ((buffer[2] != '.' && buffer[2] != ',')
- || buffer[ndigits + 2] != 'e'))
- internal_error (&dtp->common, "printf is broken"); */
+ /* Calculate total number of digits. */
+ if (ft == FMT_F)
+ ndigits = nprinted - 2;
+ else
+ ndigits = precision + 1;
/* Read the exponent back in. */
- e = atoi (&buffer[ndigits + 3]) + 1;
+ if (ft != FMT_F)
+ e = atoi (&buffer[ndigits + 3]) + 1;
+ else
+ e = 0;
/* Make sure zero comes out as 0.0e0. */
if (zero_flag)
- {
- e = 0;
- if (compile_options.sign_zero == 1)
- sign = calculate_sign (dtp, sign_bit);
- else
- sign = calculate_sign (dtp, 0);
-
- /* Handle special cases. */
- if (w == 0)
- w = d + 2;
-
- /* For this one we choose to not output a decimal point.
- F95 10.5.1.2.1 */
- if (w == 1 && ft == FMT_F)
- {
- out = write_block (dtp, w);
- if (out == NULL)
- return;
- *out = '0';
- return;
- }
-
- }
+ e = 0;
/* Normalize the fractional component. */
- buffer[2] = buffer[1];
- digits = &buffer[2];
+ if (ft != FMT_F)
+ {
+ buffer[2] = buffer[1];
+ digits = &buffer[2];
+ }
+ else
+ digits = &buffer[1];
/* Figure out where to place the decimal point. */
switch (ft)
{
case FMT_F:
- if (d == 0 && e <= 0 && dtp->u.p.scale_factor == 0)
- {
- memmove (digits + 1, digits, ndigits - 1);
- digits[0] = '0';
- e++;
+ nbefore = ndigits - precision;
+ if ((w > 0) && (nbefore > (int) size))
+ {
+ *len = w;
+ star_fill (result, w);
+ result[w] = '\0';
+ return;
}
-
- nbefore = e + dtp->u.p.scale_factor;
- if (nbefore < 0)
+ /* Make sure the decimal point is a '.'; depending on the
+ locale, this might not be the case otherwise. */
+ digits[nbefore] = '.';
+ if (p != 0)
{
- nzero = -nbefore;
- nzero_real = nzero;
- if (nzero > d)
- nzero = d;
- nafter = d - nzero;
- nbefore = 0;
+ if (p > 0)
+ {
+ memmove (digits + nbefore, digits + nbefore + 1, p);
+ digits[nbefore + p] = '.';
+ nbefore += p;
+ nafter = d;
+ nzero = 0;
+ }
+ else /* p < 0 */
+ {
+ if (nbefore + p >= 0)
+ {
+ nzero = 0;
+ memmove (digits + nbefore + p + 1, digits + nbefore + p, -p);
+ nbefore += p;
+ digits[nbefore] = '.';
+ nafter = d;
+ }
+ else
+ {
+ nzero = -(nbefore + p);
+ memmove (digits + 1, digits, nbefore);
+ nafter = d - nzero;
+ if (nafter == 0 && d > 0)
+ {
+ /* This is needed to get the correct rounding. */
+ memmove (digits + 1, digits, ndigits - 1);
+ digits[1] = '0';
+ nafter = 1;
+ nzero = d - 1;
+ }
+ else if (nafter < 0)
+ {
+ /* Reset digits to 0 in order to get correct rounding
+ towards infinity. */
+ for (i = 0; i < ndigits; i++)
+ digits[i] = '0';
+ digits[ndigits - 1] = '1';
+ nafter = d;
+ nzero = 0;
+ }
+ nbefore = 0;
+ }
+ }
}
else
{
nzero = 0;
nafter = d;
}
+
+ while (digits[0] == '0' && nbefore > 0)
+ {
+ digits++;
+ nbefore--;
+ ndigits--;
+ }
+
expchar = 0;
+ /* If we need to do rounding ourselves, get rid of the dot by
+ moving the fractional part. */
+ if (dtp->u.p.current_unit->round_status != ROUND_UNSPECIFIED
+ && dtp->u.p.current_unit->round_status != ROUND_PROCDEFINED)
+ memmove (digits + nbefore, digits + nbefore + 1, ndigits - nbefore);
break;
case FMT_E:
case FMT_D:
i = dtp->u.p.scale_factor;
- if (d <= 0 && i == 0)
+ if (d <= 0 && p == 0)
{
generate_error (&dtp->common, LIBERROR_FORMAT, "Precision not "
"greater than zero in format specifier 'E' or 'D'");
return;
}
- if (i <= -d || i >= d + 2)
+ if (p <= -d || p >= d + 2)
{
generate_error (&dtp->common, LIBERROR_FORMAT, "Scale factor "
"out of range in format specifier 'E' or 'D'");
}
if (!zero_flag)
- e -= i;
- if (i < 0)
+ e -= p;
+ if (p < 0)
{
nbefore = 0;
- nzero = -i;
- nafter = d + i;
+ nzero = -p;
+ nafter = d + p;
}
- else if (i > 0)
+ else if (p > 0)
{
- nbefore = i;
+ nbefore = p;
nzero = 0;
- nafter = (d - i) + 1;
+ nafter = (d - p) + 1;
}
- else /* i == 0 */
+ else /* p == 0 */
{
nbefore = 0;
nzero = 0;
internal_error (&dtp->common, "Unexpected format token");
}
- /* Round the value. The value being rounded is an unsigned magnitude.
- The ROUND_COMPATIBLE is rounding away from zero when there is a tie. */
+ if (zero_flag)
+ goto skip;
+
+ /* Round the value. The value being rounded is an unsigned magnitude. */
switch (dtp->u.p.current_unit->round_status)
{
+ /* For processor defined and unspecified rounding we use
+ snprintf to print the exact number of digits needed, and thus
+ let snprintf handle the rounding. On system claiming support
+ for IEEE 754, this ought to be round to nearest, ties to
+ even, corresponding to the Fortran ROUND='NEAREST'. */
+ case ROUND_PROCDEFINED:
+ case ROUND_UNSPECIFIED:
case ROUND_ZERO: /* Do nothing and truncation occurs. */
goto skip;
case ROUND_UP:
if (sign_bit)
goto skip;
- rchar = '0';
- break;
+ goto updown;
case ROUND_DOWN:
if (!sign_bit)
goto skip;
- rchar = '0';
- break;
+ goto updown;
case ROUND_NEAREST:
/* Round compatible unless there is a tie. A tie is a 5 with
all trailing zero's. */
if (digits[i] != '0')
goto do_rnd;
}
- /* It is a tie so round to even. */
+ /* It is a tie so round to even. */
switch (digits[nafter + nbefore - 1])
{
case '1':
goto skip;
}
}
- /* Fall through. */
- case ROUND_PROCDEFINED:
- case ROUND_UNSPECIFIED:
+ /* Fall through. */
+ /* The ROUND_COMPATIBLE is rounding away from zero when there is a tie. */
case ROUND_COMPATIBLE:
rchar = '5';
- /* Just fall through and do the actual rounding. */
+ goto do_rnd;
+ }
+
+ updown:
+
+ rchar = '0';
+ if (ft != FMT_F && w > 0 && d == 0 && p == 0)
+ nbefore = 1;
+ /* Scan for trailing zeros to see if we really need to round it. */
+ for(i = nbefore + nafter; i < ndigits; i++)
+ {
+ if (digits[i] != '0')
+ goto do_rnd;
}
+ goto skip;
do_rnd:
if (nbefore + nafter == 0)
+ /* Handle the case Fw.0 and value < 1.0 */
{
ndigits = 0;
- if (nzero_real == d && digits[0] >= rchar)
+ if (digits[0] >= rchar)
{
/* We rounded to zero but shouldn't have */
- nzero--;
- nafter = 1;
+ nbefore = 1;
+ digits--;
digits[0] = '1';
ndigits = 1;
}
}
else if (nbefore + nafter < ndigits)
{
- ndigits = nbefore + nafter;
- i = ndigits;
+ i = ndigits = nbefore + nafter;
if (digits[i] >= rchar)
{
/* Propagate the carry. */
skip:
/* Calculate the format of the exponent field. */
- if (expchar)
+ if (expchar && !(dtp->u.p.g0_no_blanks && e == 0))
{
edigits = 1;
for (i = abs (e); i >= 10; i /= 10)
else
edigits = 0;
- /* Zero values always output as positive, even if the value was negative
- before rounding. */
- for (i = 0; i < ndigits; i++)
+ /* Scan the digits string and count the number of zeros. If we make it
+ all the way through the loop, we know the value is zero after the
+ rounding completed above. */
+ int hasdot = 0;
+ for (i = 0; i < ndigits + hasdot; i++)
{
- if (digits[i] != '0')
+ if (digits[i] == '.')
+ hasdot = 1;
+ else if (digits[i] != '0')
break;
}
- if (i == ndigits)
+
+ /* To format properly, we need to know if the rounded result is zero and if
+ so, we set the zero_flag which may have been already set for
+ actual zero. */
+ if (i == ndigits + hasdot)
{
+ zero_flag = true;
/* The output is zero, so set the sign according to the sign bit unless
-fno-sign-zero was specified. */
if (compile_options.sign_zero == 1)
sign = calculate_sign (dtp, 0);
}
- /* Pick a field size if none was specified. */
+ /* Pick a field size if none was specified, taking into account small
+ values that may have been rounded to zero. */
if (w <= 0)
- w = nbefore + nzero + nafter + (sign != S_NONE ? 2 : 1);
-
+ {
+ if (zero_flag)
+ w = d + (sign != S_NONE ? 2 : 1) + (d == 0 ? 1 : 0);
+ else
+ {
+ w = nbefore + nzero + nafter + (sign != S_NONE ? 2 : 1);
+ w = w == 1 ? 2 : w;
+ }
+ }
+
/* Work out how much padding is needed. */
nblanks = w - (nbefore + nzero + nafter + edigits + 1);
if (sign != S_NONE)
nblanks--;
+ /* See if we have space for a zero before the decimal point. */
+ if (nbefore == 0 && nblanks > 0)
+ {
+ leadzero = 1;
+ nblanks--;
+ }
+ else
+ leadzero = 0;
+
if (dtp->u.p.g0_no_blanks)
{
w -= nblanks;
nblanks = 0;
}
- /* Create the ouput buffer. */
- out = write_block (dtp, w);
- if (out == NULL)
- return;
+ /* Create the final float string. */
+ *len = w + npad;
+ put = result;
/* Check the value fits in the specified field width. */
- if (nblanks < 0 || edigits == -1)
+ if (nblanks < 0 || edigits == -1 || w == 1 || (w == 2 && sign != S_NONE))
{
- star_fill (out, w);
+ star_fill (put, *len);
return;
}
- /* See if we have space for a zero before the decimal point. */
- if (nbefore == 0 && nblanks > 0)
- {
- leadzero = 1;
- nblanks--;
- }
- else
- leadzero = 0;
-
/* Pad to full field width. */
-
if ( ( nblanks > 0 ) && !dtp->u.p.no_leading_blank)
{
- memset (out, ' ', nblanks);
- out += nblanks;
+ memset (put, ' ', nblanks);
+ put += nblanks;
}
- /* Output the initial sign (if any). */
+ /* Set the initial sign (if any). */
if (sign == S_PLUS)
- *(out++) = '+';
+ *(put++) = '+';
else if (sign == S_MINUS)
- *(out++) = '-';
+ *(put++) = '-';
- /* Output an optional leading zero. */
+ /* Set an optional leading zero. */
if (leadzero)
- *(out++) = '0';
+ *(put++) = '0';
- /* Output the part before the decimal point, padding with zeros. */
+ /* Set the part before the decimal point, padding with zeros. */
if (nbefore > 0)
{
if (nbefore > ndigits)
{
i = ndigits;
- memcpy (out, digits, i);
+ memcpy (put, digits, i);
ndigits = 0;
while (i < nbefore)
- out[i++] = '0';
+ put[i++] = '0';
}
else
{
i = nbefore;
- memcpy (out, digits, i);
+ memcpy (put, digits, i);
ndigits -= i;
}
digits += i;
- out += nbefore;
+ put += nbefore;
}
- /* Output the decimal point. */
- *(out++) = dtp->u.p.current_unit->decimal_status == DECIMAL_POINT ? '.' : ',';
+ /* Set the decimal point. */
+ *(put++) = dtp->u.p.current_unit->decimal_status == DECIMAL_POINT ? '.' : ',';
+ if (ft == FMT_F
+ && (dtp->u.p.current_unit->round_status == ROUND_UNSPECIFIED
+ || dtp->u.p.current_unit->round_status == ROUND_PROCDEFINED))
+ digits++;
- /* Output leading zeros after the decimal point. */
+ /* Set leading zeros after the decimal point. */
if (nzero > 0)
{
for (i = 0; i < nzero; i++)
- *(out++) = '0';
+ *(put++) = '0';
}
- /* Output digits after the decimal point, padding with zeros. */
+ /* Set digits after the decimal point, padding with zeros. */
if (nafter > 0)
{
if (nafter > ndigits)
else
i = nafter;
- memcpy (out, digits, i);
+ memcpy (put, digits, i);
while (i < nafter)
- out[i++] = '0';
+ put[i++] = '0';
digits += i;
ndigits -= i;
- out += nafter;
+ put += nafter;
}
- /* Output the exponent. */
- if (expchar)
+ /* Set the exponent. */
+ if (expchar && !(dtp->u.p.g0_no_blanks && e == 0))
{
if (expchar != ' ')
{
- *(out++) = expchar;
+ *(put++) = expchar;
edigits--;
}
-#if HAVE_SNPRINTF
snprintf (buffer, size, "%+0*d", edigits, e);
-#else
- sprintf (buffer, "%+0*d", edigits, e);
-#endif
- memcpy (out, buffer, edigits);
+ memcpy (put, buffer, edigits);
+ put += edigits;
}
if (dtp->u.p.no_leading_blank)
{
- out += edigits;
- memset( out , ' ' , nblanks );
+ memset (put , ' ' , nblanks);
dtp->u.p.no_leading_blank = 0;
+ put += nblanks;
}
-#undef STR
-#undef STR1
-#undef MIN_FIELD_WIDTH
+ if (npad > 0 && !dtp->u.p.g0_no_blanks)
+ {
+ memset (put , ' ' , npad);
+ put += npad;
+ }
+
+ /* NULL terminate the string. */
+ *put = '\0';
+
+ return;
}
/* Write "Infinite" or "Nan" as appropriate for the given format. */
static void
-write_infnan (st_parameter_dt *dtp, const fnode *f, int isnan_flag, int sign_bit)
+build_infnan_string (st_parameter_dt *dtp, const fnode *f, int isnan_flag,
+ int sign_bit, char *p, size_t *len)
{
- char * p, fin;
+ char fin;
int nb = 0;
+ sign_t sign;
+ int mark;
if (f->format != FMT_B && f->format != FMT_O && f->format != FMT_Z)
{
- nb = f->u.real.w;
-
- /* If the field width is zero, the processor must select a width
- not zero. 4 is chosen to allow output of '-Inf' or '+Inf' */
-
- if (nb == 0) nb = 4;
- p = write_block (dtp, nb);
- if (p == NULL)
- return;
- if (nb < 3)
- {
- memset (p, '*',nb);
- return;
- }
+ sign = calculate_sign (dtp, sign_bit);
+ mark = (sign == S_PLUS || sign == S_MINUS) ? 8 : 7;
- memset(p, ' ', nb);
- if (!isnan_flag)
+ nb = f->u.real.w;
+ *len = nb;
+
+ /* If the field width is zero, the processor must select a width
+ not zero. 4 is chosen to allow output of '-Inf' or '+Inf' */
+
+ if ((nb == 0) || dtp->u.p.g0_no_blanks)
+ {
+ if (isnan_flag)
+ nb = 3;
+ else
+ nb = (sign == S_PLUS || sign == S_MINUS) ? 4 : 3;
+ *len = nb;
+ }
+
+ p[*len] = '\0';
+ if (nb < 3)
+ {
+ memset (p, '*', nb);
+ return;
+ }
+
+ memset(p, ' ', nb);
+
+ if (!isnan_flag)
+ {
+ if (sign_bit)
{
- if (sign_bit)
- {
-
- /* If the sign is negative and the width is 3, there is
- insufficient room to output '-Inf', so output asterisks */
-
- if (nb == 3)
- {
- memset (p, '*',nb);
- return;
- }
-
- /* The negative sign is mandatory */
-
- fin = '-';
- }
- else
-
- /* The positive sign is optional, but we output it for
- consistency */
- fin = '+';
-
- if (nb > 8)
-
- /* We have room, so output 'Infinity' */
- memcpy(p + nb - 8, "Infinity", 8);
- else
-
- /* For the case of width equals 8, there is not enough room
- for the sign and 'Infinity' so we go with 'Inf' */
- memcpy(p + nb - 3, "Inf", 3);
+ /* If the sign is negative and the width is 3, there is
+ insufficient room to output '-Inf', so output asterisks */
+ if (nb == 3)
+ {
+ memset (p, '*', nb);
+ return;
+ }
+ /* The negative sign is mandatory */
+ fin = '-';
+ }
+ else
+ /* The positive sign is optional, but we output it for
+ consistency */
+ fin = '+';
+
+ if (nb > mark)
+ /* We have room, so output 'Infinity' */
+ memcpy(p + nb - 8, "Infinity", 8);
+ else
+ /* For the case of width equals 8, there is not enough room
+ for the sign and 'Infinity' so we go with 'Inf' */
+ memcpy(p + nb - 3, "Inf", 3);
+ if (sign == S_PLUS || sign == S_MINUS)
+ {
if (nb < 9 && nb > 3)
p[nb - 4] = fin; /* Put the sign in front of Inf */
else if (nb > 8)
p[nb - 9] = fin; /* Put the sign in front of Infinity */
}
- else
- memcpy(p + nb - 3, "NaN", 3);
- return;
}
+ else
+ memcpy(p + nb - 3, "NaN", 3);
}
+}
/* Returns the value of 10**d. */
#define CALCULATE_EXP(x) \
-inline static GFC_REAL_ ## x \
+static GFC_REAL_ ## x \
calculate_exp_ ## x (int d)\
{\
int i;\
#endif
#undef CALCULATE_EXP
-/* Generate corresponding I/O format for FMT_G and output.
- The rules to translate FMT_G to FMT_E or FMT_F from DEC fortran
- LRM (table 11-2, Chapter 11, "I/O Formatting", P11-25) is:
-
- Data Magnitude Equivalent Conversion
- 0< m < 0.1-0.5*10**(-d-1) Ew.d[Ee]
- m = 0 F(w-n).(d-1), n' '
- 0.1-0.5*10**(-d-1)<= m < 1-0.5*10**(-d) F(w-n).d, n' '
- 1-0.5*10**(-d)<= m < 10-0.5*10**(-d+1) F(w-n).(d-1), n' '
- 10-0.5*10**(-d+1)<= m < 100-0.5*10**(-d+2) F(w-n).(d-2), n' '
- ................ ..........
- 10**(d-1)-0.5*10**(-1)<= m <10**d-0.5 F(w-n).0,n(' ')
- m >= 10**d-0.5 Ew.d[Ee]
-
- notes: for Gw.d , n' ' means 4 blanks
- for Gw.dEe, n' ' means e+2 blanks */
-
-#define OUTPUT_FLOAT_FMT_G(x) \
-static void \
-output_float_FMT_G_ ## x (st_parameter_dt *dtp, const fnode *f, \
- GFC_REAL_ ## x m, char *buffer, size_t size, \
- int sign_bit, bool zero_flag, int ndigits, int edigits) \
-{ \
- int e = f->u.real.e;\
- int d = f->u.real.d;\
- int w = f->u.real.w;\
- fnode *newf;\
- GFC_REAL_ ## x rexp_d;\
- int low, high, mid;\
- int ubound, lbound;\
- char *p;\
- int save_scale_factor, nb = 0;\
-\
- save_scale_factor = dtp->u.p.scale_factor;\
- newf = (fnode *) get_mem (sizeof (fnode));\
-\
- rexp_d = calculate_exp_ ## x (-d);\
- if ((m > 0.0 && m < 0.1 - 0.05 * rexp_d) || (rexp_d * (m + 0.5) >= 1.0) ||\
- ((m == 0.0) && !(compile_options.allow_std & GFC_STD_F2003)))\
- { \
- newf->format = FMT_E;\
- newf->u.real.w = w;\
- newf->u.real.d = d;\
- newf->u.real.e = e;\
- nb = 0;\
- goto finish;\
- }\
-\
- mid = 0;\
- low = 0;\
- high = d + 1;\
- lbound = 0;\
- ubound = d + 1;\
-\
- while (low <= high)\
- { \
- GFC_REAL_ ## x temp;\
- mid = (low + high) / 2;\
-\
- temp = (calculate_exp_ ## x (mid - 1) * (1 - 0.5 * rexp_d));\
-\
- if (m < temp)\
- { \
- ubound = mid;\
- if (ubound == lbound + 1)\
- break;\
- high = mid - 1;\
- }\
- else if (m > temp)\
- { \
- lbound = mid;\
- if (ubound == lbound + 1)\
- { \
- mid ++;\
- break;\
- }\
- low = mid + 1;\
- }\
- else\
- {\
- mid++;\
- break;\
- }\
- }\
-\
- if (e < 0)\
- nb = 4;\
- else\
- nb = e + 2;\
-\
- newf->format = FMT_F;\
- newf->u.real.w = f->u.real.w - nb;\
-\
- if (m == 0.0)\
- newf->u.real.d = d - 1;\
- else\
- newf->u.real.d = - (mid - d - 1);\
-\
- dtp->u.p.scale_factor = 0;\
-\
- finish:\
- output_float (dtp, newf, buffer, size, sign_bit, zero_flag, ndigits, \
- edigits);\
- dtp->u.p.scale_factor = save_scale_factor;\
-\
- free (newf);\
-\
- if (nb > 0 && !dtp->u.p.g0_no_blanks)\
- { \
- p = write_block (dtp, nb);\
- if (p == NULL)\
- return;\
- memset (p, ' ', nb);\
- }\
-}\
-
-OUTPUT_FLOAT_FMT_G(4)
-
-OUTPUT_FLOAT_FMT_G(8)
-
-#ifdef HAVE_GFC_REAL_10
-OUTPUT_FLOAT_FMT_G(10)
-#endif
-
-#ifdef HAVE_GFC_REAL_16
-OUTPUT_FLOAT_FMT_G(16)
-#endif
-
-#undef OUTPUT_FLOAT_FMT_G
-
-/* Define a macro to build code for write_float. */
+/* Define macros to build code for format_float. */
- /* Note: Before output_float is called, sprintf is used to print to buffer the
- number in the format +D.DDDDe+ddd. For an N digit exponent, this gives us
- (MIN_FIELD_WIDTH-5)-N digits after the decimal point, plus another one
- before the decimal point.
+ /* Note: Before output_float is called, snprintf is used to print to buffer the
+ number in the format +D.DDDDe+ddd.
# The result will always contain a decimal point, even if no
digits follow it
+ A sign (+ or -) always be placed before a number
- MIN_FIELD_WIDTH minimum field width
-
- * (ndigits-1) is used as the precision
+ * prec is used as the precision
e format: [-]d.ddde±dd where there is one digit before the
decimal-point character and the number of digits after it is
equal to the precision. The exponent always contains at least two
digits; if the value is zero, the exponent is 00. */
-#ifdef HAVE_SNPRINTF
-#define DTOA \
-snprintf (buffer, size, "%+-#" STR(MIN_FIELD_WIDTH) ".*" \
- "e", ndigits - 1, tmp);
+#define TOKENPASTE(x, y) TOKENPASTE2(x, y)
+#define TOKENPASTE2(x, y) x ## y
+
+#define DTOA(suff,prec,val) TOKENPASTE(DTOA2,suff)(prec,val)
+
+#define DTOA2(prec,val) \
+snprintf (buffer, size, "%+-#.*e", (prec), (val))
+
+#define DTOA2L(prec,val) \
+snprintf (buffer, size, "%+-#.*Le", (prec), (val))
+
+
+#if defined(GFC_REAL_16_IS_FLOAT128)
+#define DTOA2Q(prec,val) \
+quadmath_snprintf (buffer, size, "%+-#.*Qe", (prec), (val))
+#endif
-#define DTOAL \
-snprintf (buffer, size, "%+-#" STR(MIN_FIELD_WIDTH) ".*" \
- "Le", ndigits - 1, tmp);
+#define FDTOA(suff,prec,val) TOKENPASTE(FDTOA2,suff)(prec,val)
-#else
+/* For F format, we print to the buffer with f format. */
+#define FDTOA2(prec,val) \
+snprintf (buffer, size, "%+-#.*f", (prec), (val))
-#define DTOA \
-sprintf (buffer, "%+-#" STR(MIN_FIELD_WIDTH) ".*" \
- "e", ndigits - 1, tmp);
+#define FDTOA2L(prec,val) \
+snprintf (buffer, size, "%+-#.*Lf", (prec), (val))
-#define DTOAL \
-sprintf (buffer, "%+-#" STR(MIN_FIELD_WIDTH) ".*" \
- "Le", ndigits - 1, tmp);
+#if defined(GFC_REAL_16_IS_FLOAT128)
+#define FDTOA2Q(prec,val) \
+quadmath_snprintf (buffer, size, "%+-#.*Qf", \
+ (prec), (val))
#endif
-#define WRITE_FLOAT(x,y)\
+
+/* EN format is tricky since the number of significant digits depends
+ on the magnitude. Solve it by first printing a temporary value and
+ figure out the number of significant digits from the printed
+ exponent. Values y, 0.95*10.0**e <= y <10.0**e, are rounded to
+ 10.0**e even when the final result will not be rounded to 10.0**e.
+ For these values the exponent returned by atoi has to be decremented
+ by one. The values y in the ranges
+ (1000.0-0.5*10.0**(-d))*10.0**(3*n) <= y < 10.0*(3*(n+1))
+ (100.0-0.5*10.0**(-d))*10.0**(3*n) <= y < 10.0*(3*n+2)
+ (10.0-0.5*10.0**(-d))*10.0**(3*n) <= y < 10.0*(3*n+1)
+ are correctly rounded respectively to 1.0...0*10.0*(3*(n+1)),
+ 100.0...0*10.0*(3*n), and 10.0...0*10.0*(3*n), where 0...0
+ represents d zeroes, by the lines 279 to 297. */
+#define EN_PREC(x,y)\
{\
- GFC_REAL_ ## x tmp;\
- tmp = * (GFC_REAL_ ## x *)source;\
- sign_bit = __builtin_signbit (tmp);\
- if (!isfinite (tmp))\
- { \
- write_infnan (dtp, f, isnan (tmp), sign_bit);\
- return;\
+ volatile GFC_REAL_ ## x tmp, one = 1.0;\
+ tmp = * (GFC_REAL_ ## x *)source;\
+ if (isfinite (tmp))\
+ {\
+ nprinted = DTOA(y,0,tmp);\
+ int e = atoi (&buffer[4]);\
+ if (buffer[1] == '1')\
+ {\
+ tmp = (calculate_exp_ ## x (-e)) * tmp;\
+ tmp = one - (tmp < 0 ? -tmp : tmp);\
+ if (tmp > 0)\
+ e = e - 1;\
}\
- tmp = sign_bit ? -tmp : tmp;\
- zero_flag = (tmp == 0.0);\
-\
- DTOA ## y\
-\
- if (f->format != FMT_G)\
- output_float (dtp, f, buffer, size, sign_bit, zero_flag, ndigits, \
- edigits);\
- else \
- output_float_FMT_G_ ## x (dtp, f, tmp, buffer, size, sign_bit, \
- zero_flag, ndigits, edigits);\
+ nbefore = e%3;\
+ if (nbefore < 0)\
+ nbefore = 3 + nbefore;\
+ }\
+ else\
+ nprinted = -1;\
}\
-/* Output a real number according to its format. */
-
-static void
-write_float (st_parameter_dt *dtp, const fnode *f, const char *source, int len)
+static int
+determine_en_precision (st_parameter_dt *dtp, const fnode *f,
+ const char *source, int len)
{
+ int nprinted;
+ char buffer[10];
+ const size_t size = 10;
+ int nbefore; /* digits before decimal point - 1. */
+
+ switch (len)
+ {
+ case 4:
+ EN_PREC(4,)
+ break;
+
+ case 8:
+ EN_PREC(8,)
+ break;
-#if defined(HAVE_GFC_REAL_16) && __LDBL_DIG__ > 18
-# define MIN_FIELD_WIDTH 46
-#else
-# define MIN_FIELD_WIDTH 31
+#ifdef HAVE_GFC_REAL_10
+ case 10:
+ EN_PREC(10,L)
+ break;
#endif
-#define STR(x) STR1(x)
-#define STR1(x) #x
+#ifdef HAVE_GFC_REAL_16
+ case 16:
+# ifdef GFC_REAL_16_IS_FLOAT128
+ EN_PREC(16,Q)
+# else
+ EN_PREC(16,L)
+# endif
+ break;
+#endif
+ default:
+ internal_error (NULL, "bad real kind");
+ }
- /* This must be large enough to accurately hold any value. */
- char buffer[MIN_FIELD_WIDTH+1];
- int sign_bit, ndigits, edigits;
- bool zero_flag;
- size_t size;
+ if (nprinted == -1)
+ return -1;
- size = MIN_FIELD_WIDTH+1;
+ int prec = f->u.real.d + nbefore;
+ if (dtp->u.p.current_unit->round_status != ROUND_UNSPECIFIED
+ && dtp->u.p.current_unit->round_status != ROUND_PROCDEFINED)
+ prec += 2 * len + 4;
+ return prec;
+}
+
- /* printf pads blanks for us on the exponent so we just need it big enough
- to handle the largest number of exponent digits expected. */
- edigits=4;
+/* Generate corresponding I/O format. and output.
+ The rules to translate FMT_G to FMT_E or FMT_F from DEC fortran
+ LRM (table 11-2, Chapter 11, "I/O Formatting", P11-25) is:
- if (f->format == FMT_F || f->format == FMT_EN || f->format == FMT_G
- || ((f->format == FMT_D || f->format == FMT_E)
- && dtp->u.p.scale_factor != 0))
- {
- /* Always convert at full precision to avoid double rounding. */
- ndigits = MIN_FIELD_WIDTH - 4 - edigits;
- }
- else
- {
- /* The number of digits is known, so let printf do the rounding. */
- if (f->format == FMT_ES)
- ndigits = f->u.real.d + 1;
- else
- ndigits = f->u.real.d;
- if (ndigits > MIN_FIELD_WIDTH - 4 - edigits)
- ndigits = MIN_FIELD_WIDTH - 4 - edigits;
- }
+ Data Magnitude Equivalent Conversion
+ 0< m < 0.1-0.5*10**(-d-1) Ew.d[Ee]
+ m = 0 F(w-n).(d-1), n' '
+ 0.1-0.5*10**(-d-1)<= m < 1-0.5*10**(-d) F(w-n).d, n' '
+ 1-0.5*10**(-d)<= m < 10-0.5*10**(-d+1) F(w-n).(d-1), n' '
+ 10-0.5*10**(-d+1)<= m < 100-0.5*10**(-d+2) F(w-n).(d-2), n' '
+ ................ ..........
+ 10**(d-1)-0.5*10**(-1)<= m <10**d-0.5 F(w-n).0,n(' ')
+ m >= 10**d-0.5 Ew.d[Ee]
- switch (len)
+ notes: for Gw.d , n' ' means 4 blanks
+ for Gw.dEe, n' ' means e+2 blanks
+ for rounding modes adjustment, r, See Fortran F2008 10.7.5.2.2
+ the asm volatile is required for 32-bit x86 platforms. */
+#define FORMAT_FLOAT(x,y)\
+{\
+ int npad = 0;\
+ GFC_REAL_ ## x m;\
+ m = * (GFC_REAL_ ## x *)source;\
+ sign_bit = signbit (m);\
+ if (!isfinite (m))\
+ { \
+ build_infnan_string (dtp, f, isnan (m), sign_bit, result, res_len);\
+ return;\
+ }\
+ m = sign_bit ? -m : m;\
+ zero_flag = (m == 0.0);\
+ if (f->format == FMT_G)\
+ {\
+ int e = f->u.real.e;\
+ int d = f->u.real.d;\
+ int w = f->u.real.w;\
+ fnode newf;\
+ GFC_REAL_ ## x exp_d, r = 0.5, r_sc;\
+ int low, high, mid;\
+ int ubound, lbound;\
+ int save_scale_factor;\
+ volatile GFC_REAL_ ## x temp;\
+ save_scale_factor = dtp->u.p.scale_factor;\
+ switch (dtp->u.p.current_unit->round_status)\
+ {\
+ case ROUND_ZERO:\
+ r = sign_bit ? 1.0 : 0.0;\
+ break;\
+ case ROUND_UP:\
+ r = 1.0;\
+ break;\
+ case ROUND_DOWN:\
+ r = 0.0;\
+ break;\
+ default:\
+ break;\
+ }\
+ exp_d = calculate_exp_ ## x (d);\
+ r_sc = (1 - r / exp_d);\
+ temp = 0.1 * r_sc;\
+ if ((m > 0.0 && ((m < temp) || (r >= (exp_d - m))))\
+ || ((m == 0.0) && !(compile_options.allow_std\
+ & (GFC_STD_F2003 | GFC_STD_F2008)))\
+ || d == 0)\
+ { \
+ newf.format = FMT_E;\
+ newf.u.real.w = w;\
+ newf.u.real.d = d - comp_d;\
+ newf.u.real.e = e;\
+ npad = 0;\
+ precision = determine_precision (dtp, &newf, x);\
+ nprinted = DTOA(y,precision,m);\
+ }\
+ else \
+ {\
+ mid = 0;\
+ low = 0;\
+ high = d + 1;\
+ lbound = 0;\
+ ubound = d + 1;\
+ while (low <= high)\
+ {\
+ mid = (low + high) / 2;\
+ temp = (calculate_exp_ ## x (mid - 1) * r_sc);\
+ if (m < temp)\
+ { \
+ ubound = mid;\
+ if (ubound == lbound + 1)\
+ break;\
+ high = mid - 1;\
+ }\
+ else if (m > temp)\
+ { \
+ lbound = mid;\
+ if (ubound == lbound + 1)\
+ { \
+ mid ++;\
+ break;\
+ }\
+ low = mid + 1;\
+ }\
+ else\
+ {\
+ mid++;\
+ break;\
+ }\
+ }\
+ npad = e <= 0 ? 4 : e + 2;\
+ npad = npad >= w ? w - 1 : npad;\
+ npad = dtp->u.p.g0_no_blanks ? 0 : npad;\
+ newf.format = FMT_F;\
+ newf.u.real.w = w - npad;\
+ newf.u.real.d = m == 0.0 ? d - 1 : -(mid - d - 1) ;\
+ dtp->u.p.scale_factor = 0;\
+ precision = determine_precision (dtp, &newf, x);\
+ nprinted = FDTOA(y,precision,m);\
+ }\
+ build_float_string (dtp, &newf, buffer, size, nprinted, precision,\
+ sign_bit, zero_flag, npad, result, res_len);\
+ dtp->u.p.scale_factor = save_scale_factor;\
+ }\
+ else\
+ {\
+ if (f->format == FMT_F)\
+ nprinted = FDTOA(y,precision,m);\
+ else\
+ nprinted = DTOA(y,precision,m);\
+ build_float_string (dtp, f, buffer, size, nprinted, precision,\
+ sign_bit, zero_flag, npad, result, res_len);\
+ }\
+}\
+
+/* Output a real number according to its format. */
+
+
+static void
+get_float_string (st_parameter_dt *dtp, const fnode *f, const char *source,
+ int kind, int comp_d, char *buffer, int precision,
+ size_t size, char *result, size_t *res_len)
+{
+ int sign_bit, nprinted;
+ bool zero_flag;
+
+ switch (kind)
{
case 4:
- WRITE_FLOAT(4,)
+ FORMAT_FLOAT(4,)
break;
case 8:
- WRITE_FLOAT(8,)
+ FORMAT_FLOAT(8,)
break;
#ifdef HAVE_GFC_REAL_10
case 10:
- WRITE_FLOAT(10,L)
+ FORMAT_FLOAT(10,L)
break;
#endif
#ifdef HAVE_GFC_REAL_16
case 16:
- WRITE_FLOAT(16,L)
+# ifdef GFC_REAL_16_IS_FLOAT128
+ FORMAT_FLOAT(16,Q)
+# else
+ FORMAT_FLOAT(16,L)
+# endif
break;
#endif
default:
internal_error (NULL, "bad real kind");
}
+ return;
}