long double versions of functions in gcc under Cygwin

[lcw1964]

Greetings, all,

I am trying to port a little bit of math code to gcc, that in the
original version used the long double version of several functions (in
particular, atanl, fabsl, and expl).

I get a complie-time "unidentifi ed reference" error to the expl()
calls, but gcc seems to digest atanl and fabsl just fine. Changing expl
to exp cures the compile time problem, but I get at best double
precision in the final results. I am assuming that the use of exp() vs.
expl() is the weak link.

The GCC documentation seems to imply that expl() is supported, but I
have no idea where to find it or how to link it in properly. For that
matter, I can't seem to find prototypes in math.h for fabsl or atanl,
and they don't make gcc cough at all.

I hope this tenderfoot can find some direction, or I may resort to
singing the praises of the egregiously un-portable lcc-win32 with its
impressive 100+ digit precision qfloat library ;)

cheers,

Les

p.s. I am trying to keep this simple, so if there is a solution within
the main gcc offerings without me having to turn to the GSL, I would
like to try that first.

[Keith Thompson]

"lcw1964" <le***********@ alumni.uwo.cawr ites:

I am trying to port a little bit of math code to gcc, that in the
original version used the long double version of several functions (in
particular, atanl, fabsl, and expl).

I get a complie-time "unidentifi ed reference" error to the expl()
calls, but gcc seems to digest atanl and fabsl just fine. Changing expl
to exp cures the compile time problem, but I get at best double
precision in the final results. I am assuming that the use of exp() vs.
expl() is the weak link.

The GCC documentation seems to imply that expl() is supported, but I
have no idea where to find it or how to link it in properly. For that
matter, I can't seem to find prototypes in math.h for fabsl or atanl,
and they don't make gcc cough at all.

I hope this tenderfoot can find some direction, or I may resort to
singing the praises of the egregiously un-portable lcc-win32 with its
impressive 100+ digit precision qfloat library ;)

cheers,

Les

p.s. I am trying to keep this simple, so if there is a solution within
the main gcc offerings without me having to turn to the GSL, I would
like to try that first.

gcc is a compiler, not a complete C implementation. (Actually gcc is
a collection of compilers, but for our purposes here we can consider
only the C compiler.) The math functions are implemented by the
runtime library, not by the compiler.

In some implementations , the compiler and the runtime library are
provided together. gcc, however, generally uses whatever runtime
library is provided by the underlying operating system. On some
systems, the C runtime library happens to be one that, like gcc, is
also provided by the GNU project. I suspect you're using a system
where that isn't the case.

I suggest you ask in a newsgroup that deals with your operating system
(probably MS Windows given your mention of lcc-win32 as an
alternative).

(I don't know whether lcc-win32 provides its own C runtime library;
check the web site or ask in comp.compilers. lcc if you want more
information.)

--
Keith Thompson (The_Other_Keit h) ks***@mib.org <http://www.ghoti.net/~kst>
San Diego Supercomputer Center <* <http://users.sdsc.edu/~kst>
We must do something. This is something. Therefore, we must do this.

[lcw1964]

Keith Thompson wrote:

>
In some implementations , the compiler and the runtime library are
provided together. gcc, however, generally uses whatever runtime
library is provided by the underlying operating system. On some
systems, the C runtime library happens to be one that, like gcc, is
also provided by the GNU project. I suspect you're using a system
where that isn't the case.

I think you are right, and I was afraid it had something to do with
that.

I am learning quickly that if I wish to make the most out of higher
precision math programming, it behooves me to expand my horizons and
develop some facility with a high or even arbitrary precision package,
though I must admit that some of the source code I have contemplated
looks daunting indeed for this beginner.

In the meantime, lcc-win32 seems a reasonable, though admittedly
non-portable option. Messr. Navia's implementation of the Cephes qfloat
library seems pretty robust and at least in my barely ept hands
produces prodigious results with surprising little change to code. For
example, the very code that I am having trouble getting full long
double results with gcc/Cygwin with lcc-win32/qfloat routinely gives
102-105 digit accuracy most of the time, and at least 100 digits all of
the time. It is simply a matter of including qfloat.h, change various
commands to their qfloat equivalents (atoq, expq, atanq, etc.),
appending a q to floating point constants, and properly formatting the
output strings for printf or whatever. For my limited personal
purposes, it is about the best option I have hit upon so far, though I
do admit that cross-platform and cross-compiler compatibility would be
much more vital if my interests were less parochial.

Thanks for the feedback, though I must admit it leaves me with 3.1 gigs
of Cygwin on my hard drive that I don't know quite what to do with ;)

Les

[Dann Corbit]

"lcw1964" <le***********@ alumni.uwo.cawr ote in message
news:11******** **************@ m73g2000cwd.goo glegroups.com.. .

Greetings, all,

I am trying to port a little bit of math code to gcc, that in the
original version used the long double version of several functions (in
particular, atanl, fabsl, and expl).

I get a complie-time "unidentifi ed reference" error to the expl()
calls, but gcc seems to digest atanl and fabsl just fine. Changing expl
to exp cures the compile time problem, but I get at best double
precision in the final results. I am assuming that the use of exp() vs.
expl() is the weak link.

The GCC documentation seems to imply that expl() is supported, but I
have no idea where to find it or how to link it in properly. For that
matter, I can't seem to find prototypes in math.h for fabsl or atanl,
and they don't make gcc cough at all.

I hope this tenderfoot can find some direction, or I may resort to
singing the praises of the egregiously un-portable lcc-win32 with its
impressive 100+ digit precision qfloat library ;)

The qfloat library is by S. Moshier. You can find qfloat along with the
Cephes collection (which has tons of long double math functions) here:

http://www.moshier.net/#Cephes

cheers,

Les

p.s. I am trying to keep this simple, so if there is a solution within
the main gcc offerings without me having to turn to the GSL, I would
like to try that first.

[lcw1964]

Dann Corbit wrote:

>
The qfloat library is by S. Moshier. You can find qfloat along with the
Cephes collection (which has tons of long double math functions) here:

http://www.moshier.net/#Cephes

Thank you! I should have given Mr. Moshier proper credit. I am also
aware of the link you referred me to, my right now porting those
libraries to GCC is a little beyond my skill set, so being able to
access the qfloat functionality thru Mr. Navia's lcc-win32 "wrapper" is
a good start.

Les

[Dann Corbit]

"lcw1964" <le***********@ alumni.uwo.cawr ote in message
news:11******** **************@ n13g2000cwa.goo glegroups.com.. .

>
Dann Corbit wrote:

>>
The qfloat library is by S. Moshier. You can find qfloat along with the
Cephes collection (which has tons of long double math functions) here:

http://www.moshier.net/#Cephes

Thank you! I should have given Mr. Moshier proper credit. I am also
aware of the link you referred me to, my right now porting those
libraries to GCC is a little beyond my skill set, so being able to
access the qfloat functionality thru Mr. Navia's lcc-win32 "wrapper" is
a good start.

You don't have to know anything. They come with their own makefiles.

At most, you will have to know what kind of machine you are compiling on (if
it is not a 32 bit platform or has odd endianness or something).

The standard makefile will probably fit your situation.

Just expand this archive:
http://www.moshier.net/qlib.zip
and type "make"

The Cephes functions are even the default math functions used in some linux
distributions (IIRC).

[jaysome]

On Mon, 7 Aug 2006 21:11:08 -0700, "Dann Corbit" <dc*****@connx. com>
wrote:

>"lcw1964" <le***********@ alumni.uwo.cawr ote in message
news:11******* *************** @n13g2000cwa.go oglegroups.com. ..

>>
Dann Corbit wrote:

>>>
The qfloat library is by S. Moshier. You can find qfloat along with the
Cephes collection (which has tons of long double math functions) here:

http://www.moshier.net/#Cephes

Thank you! I should have given Mr. Moshier proper credit. I am also
aware of the link you referred me to, my right now porting those
libraries to GCC is a little beyond my skill set, so being able to
access the qfloat functionality thru Mr. Navia's lcc-win32 "wrapper" is
a good start.

You don't have to know anything. They come with their own makefiles.

At most, you will have to know what kind of machine you are compiling on (if
it is not a 32 bit platform or has odd endianness or something).

The standard makefile will probably fit your situation.

Just expand this archive:
http://www.moshier.net/qlib.zip
and type "make"

The Cephes functions are even the default math functions used in some linux
distribution s (IIRC).

There are 477 usages of "goto" in this source. That gives me a queasy
feeling.

One part of me sees me sitting through a code review and vehemently
rebuking this code after coming across about the 5th goto statement.
The other part of me sees me reviewing the black box test results that
passed and not caring about how this was coded, as long as it was
coded in Standard C. Oh the dichotomy.

--
Jay

[jacob navia]

jaysome a écrit :

There are 477 usages of "goto" in this source. That gives me a queasy
feeling.

I rewrote all the basic functions in 386 and AMD64 assembly.
The speed gain is considerable, and the gotos are even worst:

Who hasn't written a

jmp label

in assembly?

Seriously, the code is well written, and if you look at the
dates in there you will se code from eighties. And it still runs,
twenty years later.

I would like to see what code you have written in 20 years, even
if it doesn't use gotos.

Stephen Moshier has written a very good package.

One part of me sees me sitting through a code review and vehemently
rebuking this code after coming across about the 5th goto statement.

This is just dogmatic. gotos arre part of C. And they are used in
the Cephes library in a reasonable and very clear way.

The other part of me sees me reviewing the black box test results that
passed and not caring about how this was coded, as long as it was
coded in Standard C. Oh the dichotomy.

[pete]

lcw1964 wrote:

expl()

#include <float.h>

long double fs_expl(long double x);
long double fs_logl(long double x);
long double fs_sqrtl(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;

if (1 x_max) {
x_max = fs_logl(LDBL_MA X);
x_min = fs_logl(LDBL_MI N);
}
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 LDBL_EPSILON / 4);
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;

if (LDBL_MAX >= x && x 0) {
if (1 A) {
A = fs_sqrtl(2);
B = A / 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;
}

long double fs_sqrtl(long double x)
{
long int n;
long double a, b;

if (LDBL_MAX >= x && x 0) {
for (n = 0; x 2; x /= 4) {
++n;
}
while (0.5 x) {
--n;
x *= 4;
}
a = x;
b = (1 + x) / 2;
do {
x = b;
b = (a / x + x) / 2;
} while (x b);
while (n 0) {
x *= 2;
--n;
}
while (0 n) {
x /= 2;
++n;
}
} else {
if (x != 0) {
x = LDBL_MAX;
}
}
return x;
}

--
pete

[jacob navia]

pete a écrit :

lcw1964 wrote:

>>expl()

#include <float.h>

long double fs_expl(long double x);
long double fs_logl(long double x);
long double fs_sqrtl(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;

if (1 x_max) {
x_max = fs_logl(LDBL_MA X);
x_min = fs_logl(LDBL_MI N);
}
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 LDBL_EPSILON / 4);
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;

if (LDBL_MAX >= x && x 0) {
if (1 A) {
A = fs_sqrtl(2);
B = A / 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;
}

long double fs_sqrtl(long double x)
{
long int n;
long double a, b;

if (LDBL_MAX >= x && x 0) {
for (n = 0; x 2; x /= 4) {
++n;
}
while (0.5 x) {
--n;
x *= 4;
}
a = x;
b = (1 + x) / 2;
do {
x = b;
b = (a / x + x) / 2;
} while (x b);
while (n 0) {
x *= 2;
--n;
}
while (0 n) {
x /= 2;
++n;
}
} else {
if (x != 0) {
x = LDBL_MAX;
}
}
return x;
}

I find this code well DOCUMENTED isn't it?

The source of the code (who wrote it originally), the algorithms
used are well explained, the places in the code where you have to
watch for accuracy are pointed out, a nice package.

fs_sqrt is approximately 20 times slower
than the library function.