Listing 3 Computes Machine Floating-point Parameters

#include <stdio.h>
#include <math.h>
#include <float.h>

main()
{
   int beta, p;
   double a, b, eps, epsp1, sigma, nums;

   /* Discover radix */
   a = 1.0;
   do
   {
      a = 2.0 * a;
      b = a + 1.0;
   } while ((b - a) == 1.0);

   b = 1.0;
   do
      b = 2.0 * b;
   while ((a + b) == a);
   beta = (int) ((a + b) - a);
   printf("radix:\n");
   printf("\talgorithm: %d\n",beta);
   printf("\tprovided: %d\n",FLT_RADIX);

   /* Compute precision in bits */
   p = 0;
   a = 1.0;
   do
   {
      ++p;
      a *= (double) beta;
      b = a + 1.0;
   } while ((b - a) == 1.0);
   printf("precision:\n");
   printf("\talgorithm: %d\n",p);
   printf("\tprovided: %d\n",DBL_MANT_DIG);

   /* Compute machine epsilon */
   eps = 1.0;
   do
   {
      eps = 0.5 * eps;
      epsp1 = eps + 1.0;
   } while (epsp1 > 1.0);
   epsp1 = 2.0 * eps + 1.0;
   eps = epsp1 - 1.0;
   printf("machine epsilon:\n");
   printf("\talgorithm: %g\n",eps);
   printf("\tformula: %g\n",pow(FLT_RADIX,1.0-DBL_MANT_DIG));
   printf("\tprovided: %g\n",DBL_EPSILON);

   /* Compute smallest normalized magnitude */
   printf("smallest non-zero magnitude:\n");
   printf("\tformula: %g\n",pow(FLT_RADIX,DBL_MIN_EXP-1));
   printf("\tprovided: %g\n",DBL_MIN);

   /* Compute larget normalized magnitude */
   printf("largest non-zero magnitude:\n");
   printf("\tformula: %g\n",
         pow(FLT_RADIX,DBL_MAX_EXP-1) * FLT_RADIX *
         (1.0 - pow(FLT_RADIX,-DBL_MANT_DIG)));
   printf("\tprovided: %g\n",DBL_MAX);

   printf("smallest exponent: %d\n",DBL_MIN_EXP);
   printf("largest exponent: %d\n",DBL_MAX_EXP);

   nums = 2 * (FLT_RADIX - 1)
       * pow(FLT_RADIX,DBL_MANT_DIG-1)
       * (DBL_MAX_EXP - DBL_MIN_EXP + 1);
   printf("This system has %g numbers\n",nums);
   return 0;
}

/* Output (Borland C++):
radix:
   algorithm: 2
   provided: 2
precision:
   algorithm: 53
   provided: 53
machine epsilon:
   algorithm: 2.22045e-16
   formula: 2.22045e-16
   provided: 2.22045e-16
smallest non-zero magnitude:
   formula: 2.22507e-308
   provided: 2.22507e-308
largest non-zero magnitude:
   formula: 1.79769e+308
   provided: 1.79769e+308
smallest exponent: -1021
largest exponent: 1024
This system has 1.84287e+19 numbers
*/
/* End of File */