Groups | Search | Server Info | Keyboard shortcuts | Login | Register [http] [https] [nntp] [nntps]


Groups > comp.programming > #1479

Re: Exp() function reloaded

Message-ID <4F8C9952.2141@mindspring.com> (permalink)
Date 2012-04-16 18:12 -0400
From pete <pfiland@mindspring.com>
Organization PF
Newsgroups sci.math.num-analysis, alt.lang.asm, comp.programming
Subject Re: Exp() function reloaded
References <jmdbuh$fql$1@dont-email.me> <4F8A5BF2.17C5@mindspring.com> <4f8bc299$0$1377$4fafbaef@reader2.news.tin.it> <4f8bea40$0$1381$4fafbaef@reader1.news.tin.it>

Cross-posted to 3 groups.

Show all headers | View raw


io_x wrote:
> 
> "io_x" <a@b.c.invalid> ha scritto nel messaggio
> news:4f8bc299$0$1377$4fafbaef@reader2.news.tin.it...
> >
> > "pete" <pfiland@mindspring.com> ha scritto nel messaggio
> > news:4F8A5BF2.17C5@mindspring.com...
> 
> fs_expl(1000000) here would go to seg fault...
> it seems sloppy

It *is* sloppy.
I posted it only to try to make a point about the algorithm
concerning the case when there were negative arguments 
in the function call,
and I'm not sure whether or not the point applies 
to the way that hopcode is implementing the exp function.
I wrote it to be portable 
to all conforming freestanding C implementations,
as an academic exercise in C programming.

The real function is here:

http://www.mindspring.com/~pfilandr/C/fs_math/fs_math.c

#include <float.h>

long double fs_expl(long double x);
long double fs_logl(long double x);

long double fs_expl(long double x)
{
    long unsigned n, square; 
    long double b, e;
    static long double x_max, x_min, epsilon;
    static int initialized;

    if (!initialized) {
        initialized = 1;
        x_max = fs_logl(LDBL_MAX);
        x_min = fs_logl(LDBL_MIN);
        epsilon = LDBL_EPSILON / 4;
    }
    if (x_max >= x && x >= x_min) {
        for (square = 0; x > 1; x /= 2) {
            ++square;
        }
        while (-1 > x) {
            ++square;
            x /= 2;
        }
        e = b = n = 1;
        do {
            b /= n++;
            b *= x;
            e += b;
            b /= n++;
            b *= x;
            e += b;
        } while (b > epsilon);
        while (square-- != 0) {
            e *= e;
        }
    } else {
        e = x > 0 ? LDBL_MAX : 0;
    }
    return e;
}

long double fs_logl(long double x)
{
    long int n;
    long double a, b, c, epsilon;
    static long double A, B, C;
    static int initialized;
    
    if (x > 0 && LDBL_MAX >= x) {
        if (!initialized) {
            initialized = 1;
            B = 1.5;
            do {
                A = B;
                B = 1 / A + A / 2;
            } while (A > B);
            B /= 2;             
            C = fs_logl(A);
        }
        for (n = 0; x > A; x /= 2) {
            ++n;
        }
        while (B > x) {
            --n;
            x *= 2;
        }
        a = (x - 1) / (x + 1);
        x = C * n + a;
        c = a * a;
        n = 1;
        epsilon = LDBL_EPSILON * x;
        if (0 > a) {
            if (epsilon > 0) {
                epsilon = -epsilon;
            }
            do {
                n += 2;
                a *= c;
                b = a / n;
                x += b;
            } while (epsilon > b);
        } else {
            if (0 > epsilon) {
                epsilon = -epsilon;
            }
            do {
                n += 2;
                a *= c;
                b = a / n;
                x += b;
            } while (b > epsilon);
        }
        x *= 2;
    } else {
        x = -LDBL_MAX;
    }
    return x;
}


-- 
pete

Back to comp.programming | Previous | NextPrevious in thread | Find similar


Thread

Exp() function reloaded hopcode <hopcode@invalid.de> - 2012-04-15 04:31 +0200
  Re: Exp() function reloaded pete <pfiland@mindspring.com> - 2012-04-15 01:26 -0400
    Re: Exp() function reloaded hopcode <hopcode@invalid.de> - 2012-04-15 12:18 +0200
      Re: Exp() function reloaded pete <pfiland@mindspring.com> - 2012-04-15 09:30 -0400
    Re: Exp() function reloaded "io_x" <a@b.c.invalid> - 2012-04-16 09:01 +0200
      Re: Exp() function reloaded "io_x" <a@b.c.invalid> - 2012-04-16 11:50 +0200
        Re: Exp() function reloaded "io_x" <a@b.c.invalid> - 2012-04-16 15:46 +0200
          Re: Exp() function reloaded "io_x" <a@b.c.invalid> - 2012-04-16 16:07 +0200
            Re: Exp() function reloaded "io_x" <a@b.c.invalid> - 2012-04-17 09:36 +0200
              Re: Exp() function reloaded phreda <pabloreda@gmail.com> - 2012-04-19 19:43 -0700
                Re: Exp() function reloaded hopcode <hopcode@invalid.de> - 2012-04-24 14:48 +0200
                Re: Exp() function reloaded phreda <pabloreda@gmail.com> - 2012-04-25 12:36 -0700
        Re: Exp() function reloaded pete <pfiland@mindspring.com> - 2012-04-16 18:12 -0400

csiph-web