    01: #include <iostream>
    02: #include <cstdlib>
    03: #include <cmath>
    04: #include <ctime>
    05:
    06: using namespace std;
    07:
    08: /**
    09:    Sets the seed of the random number generator.
    10: */
    11: void rand_seed()
    12: {
    13:    int seed = static_cast<int>(time(0));
    14:    srand(seed);
    15: }
    16:
    17: /**
    18:    Compute a random floating point number in a range
    19:    @param a the bottom of the range
    20:    @param b the top of the range
    21:    @return a random floating point number x,
    22:    a <= x and x <= b
    23: */
    24: double rand_double(double a, double b)
    25: {
    26:    return a + (b - a) * rand() * (1.0 / RAND_MAX);
    27: }
    28:
    29: /**
    30:    Tests whether two floating-point numbers are
    31:    approximately equal.
    32:    @param x a floating-point number
    33:    @param y another floating-point number
    34:    @return true if x and y are approximately equal
    35: */
    36: bool approx_equal(double x, double y)
    37: {
    38:    const double EPSILON = 1E-14;
    39:    if (x == 0) return fabs(y) <= EPSILON;
    40:    if (y == 0) return fabs(x) <= EPSILON;
    41:    return fabs(x - y) / max(fabs(x), fabs(y)) <= EPSILON;
    42: }
    43:
    44: /* Function to be tested */
    45:
    46: /**
    47:    Computes the square root using Heron's formula
    48:    @param a an integer >= 0
    49:    @return the square root of a
    50: */
    51: double squareroot(double a)
    52: {
    53:    if (a == 0) return 0;
    54:
    55:    double xnew = a;
    56:    double xold;
    57:
    58:    do
    59:    {
    60:       xold = xnew;
    61:       xnew = (xold + a / xold) / 2;
    62:    }
    63:    while (!approx_equal(xnew, xold));
    64:
    65:    return xnew;
    66: }
    67:
    68: /* Test harness */
    69:
    70: int main()
    71: {
    72:    rand_seed();
    73:    int i;
    74:    for (i = 1; i <= 100; i++)
    75:    {
    76:       double x = rand_double(0, 1E6);
    77:       double y = squareroot(x);
    78:       if (!approx_equal(y, pow(x, 0.5)))
    79:          cout << "Test failed. ";
    80:       else
    81:          cout << "Test passed. ";
    82:       cout << "squareroot of " << x << " = " << y << "\n";
    83:    }
    84:
    85:    return 0;
    86: }