<complex.h>

[Julian V. Noble]

Dear C Mavens,

I am writing a Computing Prescription for CiSE on complex
arithmetic. I would like to be sure that I do not include
any gaffes about C99--especially about what is in the header
file <complex.h. I have Googled <complex.hand got about
177K hits, which I obviously do not wish to go through one
at a time. I also tried "<complex.hsour ce" but still got

10K hits.

Does anyone know where I can get a complete, definitive and
(more-or-less) standard listing of <complex.h? Thanks for
any URLs.
--
Julian V. Noble
Professor Emeritus of Physics
University of Virginia

[Tim Prince]

Julian V. Noble wrote:

Malcolm wrote:

>
Well, for one I would like to see how cabs and clog are implemented.

I have already found a <complex.hwhe re cabs is done naively as

cabs(z=x+iy) = sqrt(x*x+y*y)

and this is liable to overflow when, say, x is bigger in magnitude
than the square root of the largest representable double. Not too
smart. Goldberg describes how to do it correctly and I can't imagine
why someone making a library would not follow that method.

That is reasonable if it is done with promotion to the next wider
floating point format, after checking that format includes increased
exponent range, all possible at compile time. The formats were designed
specifically for such purposes. Some platforms would do it implicitly,
without any provision in source code. Evidently, then, the best
implementations are implementation dependent.
If you know the only implementation which is acceptable to you, why are
you asking for advice?

[Andrew Poelstra]

On 2006-07-09, Julian V. Noble <jv*@virginia.e duwrote:

Andrew Poelstra wrote:

>On 2006-07-09, Julian V. Noble <jv*@virginia.e duwrote:

>>Dear C Mavens,

I am writing a Computing Prescription for CiSE on complex
arithmetic. I would like to be sure that I do not include
any gaffes about C99--especially about what is in the header
file <complex.h. I have Googled <complex.hand got about
177K hits, which I obviously do not wish to go through one
at a time. I also tried "<complex.hsour ce" but still got
10K hits.
Does anyone know where I can get a complete, definitive and
(more-or-less) standard listing of <complex.h? Thanks for
any URLs.

You could try googling "C Standard Draft" to get a copy of the
latest draft (because it's free), or you could actually buy a
copy of C99. That'll give you all you need, although it will
likely be interspersed with a /lot/ of stuff you don't need.

Thanks a lot. This was very helpful. Actually I Googled C99 Standard
Draft and got WG14N1124, a draft of 2005-05-06. It is evidently
aware of the problems of overflow and branch cut location that I
cite elsewhere in this thread. Sadly, the implementer(s) of the
versions of <complex.htha t I have seen so far have not seemed
to keep the Standard's advice in mind.

It is true that there are very few standard C99 implementations out
there. I believe that the default compiler from Solaris is compliant,
as is the Coneau compiler.

There are a few others, which I can't think of at the moment.

--
Andrew Poelstra <http://www.wpsoftware. net/projects/>
To email me, use "apoelstra" at the above domain.
"You people hate mathematics." -- James Harris

[Andrew Poelstra]

On 2006-07-09, Andrew Poelstra <ap*******@loca lhost.localdoma inwrote:

On 2006-07-09, Julian V. Noble <jv*@virginia.e duwrote:

>Andrew Poelstra wrote:

>>On 2006-07-09, Julian V. Noble <jv*@virginia.e duwrote:
Dear C Mavens,

I am writing a Computing Prescription for CiSE on complex
arithmetic . I would like to be sure that I do not include
any gaffes about C99--especially about what is in the header
file <complex.h. I have Googled <complex.hand got about
177K hits, which I obviously do not wish to go through one
at a time. I also tried "<complex.hsour ce" but still got
10K hits.
Does anyone know where I can get a complete, definitive and
(more-or-less) standard listing of <complex.h? Thanks for
any URLs.
You could try googling "C Standard Draft" to get a copy of the
latest draft (because it's free), or you could actually buy a
copy of C99. That'll give you all you need, although it will
likely be interspersed with a /lot/ of stuff you don't need.

Thanks a lot. This was very helpful. Actually I Googled C99 Standard
Draft and got WG14N1124, a draft of 2005-05-06. It is evidently
aware of the problems of overflow and branch cut location that I
cite elsewhere in this thread. Sadly, the implementer(s) of the
versions of <complex.htha t I have seen so far have not seemed
to keep the Standard's advice in mind.

It is true that there are very few standard C99 implementations out
there. I believe that the default compiler from Solaris is compliant,
as is the Coneau compiler.

There are a few others, which I can't think of at the moment.

Oops! I meant "Comeau" compiler.

--
Andrew Poelstra <http://www.wpsoftware. net/projects/>
To email me, use "apoelstra" at the above domain.
"You people hate mathematics." -- James Harris

[Keith Thompson]

"Julian V. Noble" <jv*@virginia.e duwrites:

Ben Pfaff wrote:

>"Julian V. Noble" <jv*@virginia.e duwrites:

>>Does anyone know where I can get a complete, definitive and
(more-or-less) standard listing of <complex.h? Thanks for
any URLs.

Buy a copy of C99 as a PDF. It is $30 from webstore.ansi.o rg.

Thanks, but this is a one-shot, I will never have another use
for it, and no one is going to reimburse me.

N1124 is probably better for your purposes anyway. It includes all of
C99 with two Technical Corrigenda merged in and marked with change
bars.

--
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.

[Julian V. Noble]

Tim Prince wrote:

Julian V. Noble wrote:

>Malcolm wrote:

>>
Well, for one I would like to see how cabs and clog are implemented.

I have already found a <complex.hwhe re cabs is done naively as

cabs(z=x+iy) = sqrt(x*x+y*y)

and this is liable to overflow when, say, x is bigger in magnitude
than the square root of the largest representable double. Not too
smart. Goldberg describes how to do it correctly and I can't imagine
why someone making a library would not follow that method.

That is reasonable if it is done with promotion to the next wider
floating point format, after checking that format includes increased
exponent range, all possible at compile time. The formats were designed
specifically for such purposes. Some platforms would do it implicitly,
without any provision in source code. Evidently, then, the best
implementations are implementation dependent.
If you know the only implementation which is acceptable to you, why are
you asking for advice?

Wasn't asking for advice. Was asking for source for #include file. So
I could explain pitfalls of complex libraries to people who might not
be aware of them. Was also not singling out C, I will be looking also
at Python.

--
Julian V. Noble
Professor Emeritus of Physics
University of Virginia

[Julian V. Noble]

Rouben Rostamian wrote:

In article <e8**********@m urdoch.acc.Virg inia.EDU>,
Julian V. Noble <jv*@virginia.e duwrote:

>Well, for one I would like to see how cabs and clog are implemented.

I have already found a <complex.hwhe re cabs is done naively as

cabs(z=x+iy) = sqrt(x*x+y*y)

Section 7.3.4 of the C99 Standard addresses exactly this issue.
It implies, and I paraphrase, that the implementation may have
its own clever way of computing cabs, however you can tell it
to use the formula the you have cited using a preprocessor
switch as in:

#include <complex.h>
#pragma STDC CX_LIMITED_RANG E on_or_off_switc h

The default value of the switch is off. Chances are that
the formula that you found in complex.h is guarded by a
conditional test of CX_LIMITED_RANG E.

>The complex log and related functions have implicit decisions about
where the branch line is located.

The C99 Standard is explicit about the location of branch cuts.
For instance for clog it says the branch cut is along the
negative real axis. An implementation which puts the cut anywhere
else is not standard-compliant.

Yes it is explicit. I found a draft version eventually after a kind
soul here steered me to it. The point is, one must _know_ where the
library puts the branch cuts. It is not inherently obvious to someone
who just writes a program. And different languages may make different
choices--AFIK there is no overall standard for this equivalent to
IEEE 754 for binary floating point operations.

For further details, see the C99 Standard. You may find a printed
copy in your library. See:

http://tinyurl.com/ns9o5

Thanks for the help.

--
Julian V. Noble
Professor Emeritus of Physics
University of Virginia

[Keith Thompson]

"Julian V. Noble" <jv*@virginia.e duwrites:
[...]

Wasn't asking for advice. Was asking for source for #include file. So
I could explain pitfalls of complex libraries to people who might not
be aware of them. Was also not singling out C, I will be looking also
at Python.

Are you looking for the standard's definition of what's required to be
in <complex.h>, or are you looking for the actual contents of the
header?

A standard header for a given implementation can, and often does,
contain all sorts of implementation-specific stuff. Usually a
programmer shouldn't care about the details. (An example of this is
the definition of type FILE in <stdio.h>.)

The description in the standard, on the other hand, establishes a
contract between the implementer and the programmer.

(A minor point: the standard headers aren't necessarily implemented as
files; they can be built into the compiler or stored in some other
form.)

--
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.

[P.J. Plauger]

"Julian V. Noble" <jv*@virginia.e duwrote in message
news:e8******** **@murdoch.acc. Virginia.EDU...

Wasn't asking for advice. Was asking for source for #include file. So
I could explain pitfalls of complex libraries to people who might not
be aware of them. Was also not singling out C, I will be looking also
at Python.

Well, here's some free advice anyway. You won't find the worst problems
presented in the header file. Yes, some incredibly naive implementations
will define cabs as the inline you describe, but the real issues with
complex arithmetic are:

-- avoiding intermediate overflows

-- computing accurate results in sensitive regions, such as near
zeros

-- getting all the Inf and NaN cases right

The branch cuts tend to take care of themselves, as most implementers
use the same old approximations that Fortran has relied on for half a
century.

Now that our new product line has been released, we'll be adding our
comparative analysis of C libraries over the coming weeks. You'll find
that complex math libraries have some surprising soft spots (like clog
near the unit circle) which tend to be of more practical importance
than bogus overflows on cabs calls.

P.J. Plauger
Dinkumware, Ltd.
http://www.dinkumware.com