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


Groups > comp.programming > #1479

Re: Exp() function reloaded

Path csiph.com!newsfeed.hal-mli.net!feeder3.hal-mli.net!newsfeed.hal-mli.net!feeder1.hal-mli.net!border3.nntp.dca.giganews.com!Xl.tags.giganews.com!border1.nntp.dca.giganews.com!nntp.giganews.com!local2.nntp.dca.giganews.com!nntp.earthlink.com!news.earthlink.com.POSTED!not-for-mail
NNTP-Posting-Date Mon, 16 Apr 2012 17:12:32 -0500
Message-ID <4F8C9952.2141@mindspring.com> (permalink)
Date Mon, 16 Apr 2012 18:12:34 -0400
From pete <pfiland@mindspring.com>
Reply-To pfiland@mindspring.com
Organization PF
X-Mailer Mozilla 3.04Gold (WinNT; I)
MIME-Version 1.0
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>
Content-Type text/plain; charset=us-ascii
Content-Transfer-Encoding 7bit
Lines 130
X-Usenet-Provider http://www.giganews.com
NNTP-Posting-Host 4.154.221.186
X-Trace sv3-FtV4XxG/+0tu78uxnayKsPMKXzs1ay6mKw0IXa3tjO9+ciS7sWTJxtF6/IdbDyAzCsm7fxEZ6IyOgME!NhVz0U8mXSojiV3nonpJZT4TpuBgGTBNSCVe0tO/HmEIojf4CHfSDD+dP9o7VHPb3b5+Pm/aG2rf!MDbmEX8qQQU=
X-Abuse-and-DMCA-Info Please be sure to forward a copy of ALL headers
X-Abuse-and-DMCA-Info Otherwise we will be unable to process your complaint properly
X-Postfilter 1.3.40
X-Original-Bytes 4296
Xref csiph.com comp.programming:1479

Cross-posted to 3 groups.

Show key headers only | 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