...

> > I have confidence in the square root finding while (B > A) loop,

> > but the while (x != b), has me concerned.

> > double ln(double x)

> > {

> > int n;

> > double a, b, c;

> > static double A, B, C;

> > if (!A) {

> So in effect you're depending upon all static variables being

> initialized to zero, or in this case, 0.0?

6.7.8p10: "If an object that has static storage duration is not

initialized explicitly, then: ... if it has arithmetic type, it is

initialized to (positive or unsigned) zero;"

Quote:> > A = 1.5;

> > do {

> > B = A;

> > A = 1 / B + B / 2;

> > } while (B > A); /* A = sqrt(2) */

> > B /= 2; /* B = sqrt(2) / 2 */

> > C = ln(A); /* C = ln(sqrt(2)) */

> And here you're initializing C with a recursive call to the remainder of

> the same function, which is is going to use C, except that it is

> currently 0.0? Does that even work? :-)

Of course it does. The call to ln() is recursive only if A is 0.0. By

the time the recursive call is made, A has already been set to a

non-zero value, so the recursion terminates with that second call. C

will have a value of 0.0 during that recursive call, which apparantly

doesn't cause a problem. It acquires the correct value of ln(sqrt(2))

as a result of evaluating the call.

I can't vouch for the validity of the algorithm being used (actually,

I could, but I don't want to spend the time), but it does indeed give

correct results for the values I tried (1.0, 2.0, and 4.0).

Quote:> > }

> > for (n = 0; x > A; x /= 2) {

> > ++n;

> > }

> > while (B > x) {

> > --n;

> > x *= 2;

> > }

> > a = (x - 1) / (x + 1);

> > b = C * n + a;

> So on the first (recursive) call, when you're trying to determine

> ln(sqrt(2)), C will be 0.0, right?

> What's wrong with simply stating this directly instead?

> static double A = 1.4142135623730951,

> B = 0.70710678118654757,

> C = 0.3465735902799727;

Nothing. These are two different valid ways to do it. The method he

gave only produces the wrong result if the code itself is wrong; your

approach also gives a wrong result if the initialization of the static

data is wrong; it's debateable which form is more reliable. Yours is

certainly faster, but is less portable, because you'll need to add

additional digits when porting to machines with smaller values for

DBL_EPSILON. For C99 code, hexadecimal constants would be preferred;

this is precisely the kind of situation for which hexadecimal floating

point constants were invented.