Sine code for ANSI C

[P.J. Plauger]

Re: Sine code for ANSI C

"Dan Pop" <Dan.Pop@cern.c h> wrote in message
news:c8fkhr$7u4 $9@sunnews.cern .ch...
[color=blue]
> In <V2tqc.138401$G _.102626@nwrddc 02.gnilink.net> "P.J. Plauger"[/color]
<pjp@dinkumware .com> writes:[color=blue]
>[color=green]
> >"Dan Pop" <Dan.Pop@cern.c h> wrote in message
> >news:c8dd2g$6d g$15@sunnews.ce rn.ch...
> >[color=darkred]
> >> In <40A8F445.DBE0F 395@yahoo.com> CBFalconer <cbfalconer@yah oo.com>[/color][/color][/color]
writes:[color=blue][color=green][color=darkred]
> >>
> >> >Dan Pop wrote:
> >> >>
> >> >... snip ...
> >> >>
> >> >> The point is not whether it can be calculated, but rather how much
> >> >> precision should the calculation produce? Does it makes sense to
> >> >> compute sin(DBL_MAX) with 53-bit precision, ignoring the fact that
> >> >> DBL_MAX stands for an interval so large as to make this function
> >> >> call completely devoid of any meaning?
> >> >
> >> >PJPs point is that that interval is the range of values specified,
> >> >with the actual value being the most likely value in the range.
> >> >The shape of the probability curve in the interval is known only
> >> >to the end user. It may be an impulse function precisely at the
> >> >specified value.
> >>
> >> Or it may not. If I deserve 5 digits of precision, why should the
> >> implementation waste precious CPU resources to give me 16 digits of
> >> *bogus* precision?[/color]
> >
> >Because it doesn't know that it shouldn't? You could supply a useful
> >hint by calling sinf, in this particular case at least.[/color]
>
> sinf is useless in portable C code.[/color]

But wait a minute. You've just outlined a "concrete example" where
you know you deserve only 5 digits of precision. sinf will meet
your needs *portably*. And it won't waste your precious CPU resources
to give you 16 digits of *bogus* precison. How can you call that
useless when you just outlined a concrete use case where it portably
meets your needs exactly?
[color=blue][color=green][color=darkred]
> >> To deserve 16 digits of precision, it is *my* job to call sin() with
> >> an argument of the right magnitude. It's as simple as that.[/color]
> >
> >Simple for you, perhaps, but the interface for sin currently fails
> >to provide any information about how many of the bits you're supplying
> >should be taken seriously.[/color]
>
> ALL of them, as I have already explained you. It's the magnitude of the
> value that decides how many digits of precision I deserve in the result.[/color]

But you've yet to commit to any *concrete* rule for deciding how many
digits you deserve given the magnitude of the argument value. And you
say we're at fault if a library takes more time than it should. We're
looking to you for guidance in this area, particularly since you seem
to have such a clear vision of the right way to do things. Please provide
wording for a correction to the C Standard so we can all critique it and
then adopt it. We're counting on you for leadership.

P.J. Plauger
Dinkumware, Ltd.

http://www.dinkumware.com

[Dan Pop]

Re: Sine code for ANSI C

In <Hxysr2.143@cwi .nl> "Dik T. Winter" <Dik.Winter@cwi .nl> writes:
[color=blue]
>In article <c8fl7p$7u4$11@ sunnews.cern.ch > Dan.Pop@cern.ch (Dan Pop) writes:[color=green]
> > In <HxyHxo.Hxo@cwi .nl> "Dik T. Winter" <Dik.Winter@cwi .nl> writes:[/color]
>...[color=green][color=darkred]
> > > > I would certainly look carefully at them, especially
> > > > if done with matrix inversion.
> > >
> > >Also with Gauss' elimination. If you really think that the input
> > >numbers are intervals, in almost all cases the result makes no sense.[/color]
> >
> > Unless the input numbers can be exactly represented in floating point,
> > it is sheer foolishness to consider them as exact values. Furthermore,
> > in real life applications, most input numbers are obtained with a
> > precision even lower than that provided by the floating point
> > representation, so they stand for even larger intervals. That's why
> > the results are intervals, too. Consider the difference between
> > 1 +/- 1000 and 1 +/- 0.001. Only a patent fool will display the result
> > as a plain 1 in both cases.[/color]
>
>Ah, so, what method do *you* use to solve linear systems on the computer?[/color]

I don't solve linear systems on the computer. Yet, as far as I know, all
the usual methods come bundled with an error analysis, that allows
computing the precision of the results. If you don't compute this
precision, what's the point in solving the linear systems in the first
place? There are more efficient random number generators...

Dan
--
Dan Pop
DESY Zeuthen, RZ group
Email: Dan.Pop@ifh.de

[Dan Pop]

Re: Sine code for ANSI C

In <EJKqc.512$r6.1 42@nwrddc03.gni link.net> "P.J. Plauger" <pjp@dinkumware .com> writes:
[color=blue]
>"Dan Pop" <Dan.Pop@cern.c h> wrote in message
>news:c8fkhr$7u 4$9@sunnews.cer n.ch...
>[color=green]
>> In <V2tqc.138401$G _.102626@nwrddc 02.gnilink.net> "P.J. Plauger"[/color]
><pjp@dinkumwar e.com> writes:[color=green]
>>[color=darkred]
>> >"Dan Pop" <Dan.Pop@cern.c h> wrote in message
>> >news:c8dd2g$6d g$15@sunnews.ce rn.ch...
>> >
>> >> In <40A8F445.DBE0F 395@yahoo.com> CBFalconer <cbfalconer@yah oo.com>[/color][/color]
>writes:[color=green][color=darkred]
>> >>
>> >> >Dan Pop wrote:
>> >> >>
>> >> >... snip ...
>> >> >>
>> >> >> The point is not whether it can be calculated, but rather how much
>> >> >> precision should the calculation produce? Does it makes sense to
>> >> >> compute sin(DBL_MAX) with 53-bit precision, ignoring the fact that
>> >> >> DBL_MAX stands for an interval so large as to make this function
>> >> >> call completely devoid of any meaning?
>> >> >
>> >> >PJPs point is that that interval is the range of values specified,
>> >> >with the actual value being the most likely value in the range.
>> >> >The shape of the probability curve in the interval is known only
>> >> >to the end user. It may be an impulse function precisely at the
>> >> >specified value.
>> >>
>> >> Or it may not. If I deserve 5 digits of precision, why should the
>> >> implementation waste precious CPU resources to give me 16 digits of
>> >> *bogus* precision?
>> >
>> >Because it doesn't know that it shouldn't? You could supply a useful
>> >hint by calling sinf, in this particular case at least.[/color]
>>
>> sinf is useless in portable C code.[/color]
>
>But wait a minute. You've just outlined a "concrete example" where
>you know you deserve only 5 digits of precision. sinf will meet
>your needs *portably*. And it won't waste your precious CPU resources
>to give you 16 digits of *bogus* precison. How can you call that
>useless when you just outlined a concrete use case where it portably
>meets your needs exactly?[/color]

sinf is a C99 invention and the portability of C99 code is extremely
impaired by the number of implementations claiming C99-conformance.

C99 is perfectly futile to anyone concerned about the portability of his
code. At best, it's a guide about what C89 user space identifiers not
to use in your code (because C99 has no qualms stomping on the C89 user
name space).
[color=blue][color=green][color=darkred]
>> >> To deserve 16 digits of precision, it is *my* job to call sin() with
>> >> an argument of the right magnitude. It's as simple as that.
>> >
>> >Simple for you, perhaps, but the interface for sin currently fails
>> >to provide any information about how many of the bits you're supplying
>> >should be taken seriously.[/color]
>>
>> ALL of them, as I have already explained you. It's the magnitude of the
>> value that decides how many digits of precision I deserve in the result.[/color]
>
>But you've yet to commit to any *concrete* rule for deciding how many
>digits you deserve given the magnitude of the argument value.[/color]

I have already provided the rule, last week. And when I asked if it was
concrete enough for you, you didn't object.

Dan
--
Dan Pop
DESY Zeuthen, RZ group
Email: Dan.Pop@ifh.de

[P.J. Plauger]

Re: Sine code for ANSI C

"Dan Pop" <Dan.Pop@cern.c h> wrote in message
news:c8fu88$im7 $1@sunnews.cern .ch...
[color=blue][color=green][color=darkred]
> >> sinf is useless in portable C code.[/color]
> >
> >But wait a minute. You've just outlined a "concrete example" where
> >you know you deserve only 5 digits of precision. sinf will meet
> >your needs *portably*. And it won't waste your precious CPU resources
> >to give you 16 digits of *bogus* precison. How can you call that
> >useless when you just outlined a concrete use case where it portably
> >meets your needs exactly?[/color]
>
> sinf is a C99 invention and the portability of C99 code is extremely
> impaired by the number of implementations claiming C99-conformance.[/color]

Actually it's described in C89.
[color=blue]
> C99 is perfectly futile to anyone concerned about the portability of his
> code. At best, it's a guide about what C89 user space identifiers not
> to use in your code (because C99 has no qualms stomping on the C89 user
> name space).[/color]

sinf has been reserved since 1989, so I can't imagine how using it
could possibly "stomp on" any program. But that's okay, I understand
that you're trolling.
[color=blue][color=green][color=darkred]
> >> >> To deserve 16 digits of precision, it is *my* job to call sin() with
> >> >> an argument of the right magnitude. It's as simple as that.
> >> >
> >> >Simple for you, perhaps, but the interface for sin currently fails
> >> >to provide any information about how many of the bits you're supplying
> >> >should be taken seriously.
> >>
> >> ALL of them, as I have already explained you. It's the magnitude of[/color][/color][/color]
the[color=blue][color=green][color=darkred]
> >> value that decides how many digits of precision I deserve in the[/color][/color][/color]
result.[color=blue][color=green]
> >
> >But you've yet to commit to any *concrete* rule for deciding how many
> >digits you deserve given the magnitude of the argument value.[/color]
>
> I have already provided the rule, last week. And when I asked if it was
> concrete enough for you, you didn't object.[/color]

Oh sorry, I thought you were explaining your (mis)understand ing
of the problem with this interaction:

: > >So all the sine function needs is an indication of how much fuzz to
: > >assume in the argument. I've yet to hear a specific criterion. Concrete
: > >example please.
: >
: > Assume all the bits of the argument are exact and the *only* source
: > of fuziness is the magnitude of the value. By the time the least
: > significant bit of the mantissa means more than 2 * pi, the fuziness
: > is complete (any result in the -1..1 range is correct).
: >
: > Is that concrete enough for you?
:
: Yes, it's concrete enough to show that you still don't get it.

So what you're committing to is the requirement that the sine
should *not* be computed for any value where 1 ulp weighs more
than 2*pi. And presumably, sine should deliver the best possible
approximation to any argument less than this magnitude. (This
is fine with me, since it still requires an argument reduction
using about 2 1/2 words of precision, so it will break the
vast majority of sine functions that have been in use for the
past several decades. But hey, we're covered, so I'm happy.)
Of course, you still need to specify what to return for the
garbage value, and what error to report (if any). I await your
further guidance.

So you're willing to stand behind this criterion even after I
described, in a later posting, the various ways in which it was
arbitrary? If this is submitted as a DR to WG14, you're
willing, as it's original author and proposer, to defend it
against cranky attacks by people who haven't thought through
the matter as carefully as you?

Brave of you. We can call it the Pop Criterion.

P.J. Plauger
Dinkumware, Ltd.

http://www.dinkumware.com

[Christian Bau]

Re: Sine code for ANSI C

In article <c8fsdo$ces$3@s unnews.cern.ch> , Dan.Pop@cern.ch (Dan Pop)
wrote:
[color=blue]
> I don't solve linear systems on the computer. Yet, as far as I know, all
> the usual methods come bundled with an error analysis, that allows
> computing the precision of the results. If you don't compute this
> precision, what's the point in solving the linear systems in the first
> place? There are more efficient random number generators...[/color]

It depends. Do you want a set of numbers that is close to the correct
solution of the systems, or do you want a set of numbers that is very
close to being a solution, that is each equation is solved with a small
error?

For many systems, you can get the second much easier than the first. Or
lets say you want to minimise f (x). It often doesn't matter whether you
find an x' that is close to the x with the smallest value or not, but
what you want is an x' such that f (x') is close to the smallest value
of f (x).

[Dik T. Winter]

Re: Sine code for ANSI C

In article <christian.ba u-759E4A.20583019 052004@slb-newsm1.svr.pol. co.uk> Christian Bau <christian.bau@ cbau.freeserve. co.uk> writes:[color=blue]
> In article <c8fsdo$ces$3@s unnews.cern.ch> , Dan.Pop@cern.ch (Dan Pop)
> wrote:[color=green]
> > I don't solve linear systems on the computer. Yet, as far as I know, all
> > the usual methods come bundled with an error analysis, that allows
> > computing the precision of the results. If you don't compute this
> > precision, what's the point in solving the linear systems in the first
> > place? There are more efficient random number generators...[/color]
>
> It depends. Do you want a set of numbers that is close to the correct
> solution of the systems, or do you want a set of numbers that is very
> close to being a solution, that is each equation is solved with a small
> error?[/color]

The actual methods most commonly used give the exact solution to a system
that differs from the original system with a small error. This does *not*
mean that the solution is close to the exact solution of the exact original
problem. It is close to it when the condition number of the original
matrix is good. And that is what the error analysis does do (backward
error analysis as introduced by J.H. Wilkinson).

To compute the precision of the result you have to do much more than is
feasable in most cases. Have a look at Karlsruhe arithmetic from the
University at Karlsruhe.
--
dik t. winter, cwi, kruislaan 413, 1098 sj amsterdam, nederland, +31205924131
home: bovenover 215, 1025 jn amsterdam, nederland; http://www.cwi.nl/~dik/

[Dan Pop]

Re: Sine code for ANSI C

In <HxzJFy.HxG@cwi .nl> "Dik T. Winter" <Dik.Winter@cwi .nl> writes:
[color=blue]
>To compute the precision of the result you have to do much more than is
>feasable in most cases.[/color]

What is the value of a result with unknown precision? To repeat my
earlier example, if you don't know if it's 1 +/- 0.001 or 1 +/- 1000,
the value of 1 is of pretty little use.

Dan
--
Dan Pop
DESY Zeuthen, RZ group
Email: Dan.Pop@ifh.de

[Dan Pop]

Re: Sine code for ANSI C

In <9LMqc.1021$r6. 1019@nwrddc03.g nilink.net> "P.J. Plauger" <pjp@dinkumware .com> writes:
[color=blue]
>"Dan Pop" <Dan.Pop@cern.c h> wrote in message
>news:c8fu88$im 7$1@sunnews.cer n.ch...
>[color=green][color=darkred]
>> >> sinf is useless in portable C code.
>> >
>> >But wait a minute. You've just outlined a "concrete example" where
>> >you know you deserve only 5 digits of precision. sinf will meet
>> >your needs *portably*. And it won't waste your precious CPU resources
>> >to give you 16 digits of *bogus* precison. How can you call that
>> >useless when you just outlined a concrete use case where it portably
>> >meets your needs exactly?[/color]
>>
>> sinf is a C99 invention and the portability of C99 code is extremely
>> impaired by the number of implementations claiming C99-conformance.[/color]
>
>Actually it's described in C89.[/color]

Are you trolling? It may be described in C89 (as a reserved identifier)
but no conforming C89 implementation is required to provide it. Your
own published implementation of the standard C89 library doesn't...
[color=blue][color=green]
>> C99 is perfectly futile to anyone concerned about the portability of his
>> code. At best, it's a guide about what C89 user space identifiers not
>> to use in your code (because C99 has no qualms stomping on the C89 user
>> name space).[/color]
>
>sinf has been reserved since 1989, so I can't imagine how using it
>could possibly "stomp on" any program. But that's okay, I understand
>that you're trolling.[/color]

You're trolling. I was talking in general terms: there are dozens of C99
standard library functions whose names belong to the C89 user name space:
all of <complex.h> and <fenv.h>, as well as additions to other headers,
including something as innocuous as round().
[color=blue][color=green][color=darkred]
>> >> >> To deserve 16 digits of precision, it is *my* job to call sin() with
>> >> >> an argument of the right magnitude. It's as simple as that.
>> >> >
>> >> >Simple for you, perhaps, but the interface for sin currently fails
>> >> >to provide any information about how many of the bits you're supplying
>> >> >should be taken seriously.
>> >>
>> >> ALL of them, as I have already explained you. It's the magnitude of[/color][/color]
>the[color=green][color=darkred]
>> >> value that decides how many digits of precision I deserve in the[/color][/color]
>result.[color=green][color=darkred]
>> >
>> >But you've yet to commit to any *concrete* rule for deciding how many
>> >digits you deserve given the magnitude of the argument value.[/color]
>>
>> I have already provided the rule, last week. And when I asked if it was
>> concrete enough for you, you didn't object.[/color]
>
>Oh sorry, I thought you were explaining your (mis)understand ing
>of the problem with this interaction:
>
>: > >So all the sine function needs is an indication of how much fuzz to
>: > >assume in the argument. I've yet to hear a specific criterion. Concrete
>: > >example please.
>: >
>: > Assume all the bits of the argument are exact and the *only* source
>: > of fuziness is the magnitude of the value. By the time the least
>: > significant bit of the mantissa means more than 2 * pi, the fuziness
>: > is complete (any result in the -1..1 range is correct).
>: >
>: > Is that concrete enough for you?
>:
>: Yes, it's concrete enough to show that you still don't get it.
>
>So what you're committing to is the requirement that the sine
>should *not* be computed for any value where 1 ulp weighs more
>than 2*pi.[/color]

That's too strong. I would still expect

sin(DBL_MAX) * sin(DBL_MAX) + cos(DBL_MAX) * cos(DBL_MAX)

to be reasonably close to 1.0, if not equal to it. But, other than that,
any value in the [-1, 1] range is OK once 1 ulp exceeds 2*pi.
[color=blue]
>And presumably, sine should deliver the best possible
>approximatio n to any argument less than this magnitude.[/color]

Not even that. If the value in question covers the range [min, max],
sine should provide the best possible approximation of sin(x), where x
is one value in the range [min, max].
[color=blue]
>(This
>is fine with me, since it still requires an argument reduction
>using about 2 1/2 words of precision, so it will break the
>vast majority of sine functions that have been in use for the
>past several decades.[/color]

Only the ones that have been *misused*. If I measure the volume of a
sphere with a precision of 1 cubical centimeter, does it make any sense
to compute *and display* its radius with picometer precision?
[color=blue]
>But hey, we're covered, so I'm happy.)
>Of course, you still need to specify what to return for the
>garbage value, and what error to report (if any). I await your
>further guidance.[/color]

For the garbage value, see above. Set errno to any value you consider
sensible and document it. ERANGE would be fine with me.
[color=blue]
>So you're willing to stand behind this criterion even after I
>described, in a later posting, the various ways in which it was
>arbitrary?[/color]

Yup. I didn't buy any of your sophistry.
[color=blue]
>If this is submitted as a DR to WG14, you're
>willing, as it's original author and proposer, to defend it
>against cranky attacks by people who haven't thought through
>the matter as carefully as you?[/color]

Considering that the current standard doesn't impose *any* precision on
sin, I would not expect my DR to be taken seriously. Whatever
considerations led to the current committee decision would be still valid
after the submission of my DR (no precision specification is *more*
permissive than *any* precision specification).

If I were to submit a DR, it would be one requiring each implementation
to document the precision of the 4 arithmetic operations on each floating
point type. But I know, from past c.s.c discussion, that the committee
would turn a deaf ear to my DR.
[color=blue]
>Brave of you. We can call it the Pop Criterion.[/color]

Call it the common sense criterion.

Dan
--
Dan Pop
DESY Zeuthen, RZ group
Email: Dan.Pop@ifh.de

[P.J. Plauger]

Re: Sine code for ANSI C

"Dan Pop" <Dan.Pop@cern.c h> wrote in message
news:c8kvub$34v $17@sunnews.cer n.ch...
[color=blue][color=green]
> >: > Assume all the bits of the argument are exact and the *only* source
> >: > of fuziness is the magnitude of the value. By the time the least
> >: > significant bit of the mantissa means more than 2 * pi, the fuziness
> >: > is complete (any result in the -1..1 range is correct).
> >: >
> >: > Is that concrete enough for you?
> >:
> >: Yes, it's concrete enough to show that you still don't get it.
> >
> >So what you're committing to is the requirement that the sine
> >should *not* be computed for any value where 1 ulp weighs more
> >than 2*pi.[/color]
>
> That's too strong. I would still expect
>
> sin(DBL_MAX) * sin(DBL_MAX) + cos(DBL_MAX) * cos(DBL_MAX)
>
> to be reasonably close to 1.0, if not equal to it. But, other than that,
> any value in the [-1, 1] range is OK once 1 ulp exceeds 2*pi.[/color]

There you go again, getting un-concrete. You'll have to specify
what "reasonably " means, or some nit-picker will take apart
your revision of the sine specification.
[color=blue][color=green]
> >And presumably, sine should deliver the best possible
> >approximatio n to any argument less than this magnitude.[/color]
>
> Not even that. If the value in question covers the range [min, max],
> sine should provide the best possible approximation of sin(x), where x
> is one value in the range [min, max].[/color]

Okay, then I see no reason why we shouldn't apply this principle
to *every* math function. That would mean, for example, that
the range of permissible values can be huge for large arguments
to exp, lgamma, tgamma, as well as many combinations of arguments
to pow.
[color=blue][color=green]
> >(This
> >is fine with me, since it still requires an argument reduction
> >using about 2 1/2 words of precision, so it will break the
> >vast majority of sine functions that have been in use for the
> >past several decades.[/color]
>
> Only the ones that have been *misused*. If I measure the volume of a
> sphere with a precision of 1 cubical centimeter, does it make any sense
> to compute *and display* its radius with picometer precision?[/color]

No. But if you specify the behavior of a math function and don't
provide a way to specify the uncertainty in its arguments, does
it make any sense to demand that the math functions know how
much uncertainty is there? I've been talking about the requirements
on the authors of math functions and you keep going back to how
they should be properly used. Two different universes of discourse.
[color=blue][color=green]
> >But hey, we're covered, so I'm happy.)
> >Of course, you still need to specify what to return for the
> >garbage value, and what error to report (if any). I await your
> >further guidance.[/color]
>
> For the garbage value, see above. Set errno to any value you consider
> sensible and document it. ERANGE would be fine with me.[/color]

Okay, so a programmer knows that a call to sine was silly if

1) the return value is between -1 and +1, and

2) errno indicates that some function value is either too large to
represent (overflow) or too small to represent (underflow).

That should stand up well to criticism.
[color=blue][color=green]
> >So you're willing to stand behind this criterion even after I
> >described, in a later posting, the various ways in which it was
> >arbitrary?[/color]
>
> Yup. I didn't buy any of your sophistry.
>[color=green]
> >If this is submitted as a DR to WG14, you're
> >willing, as it's original author and proposer, to defend it
> >against cranky attacks by people who haven't thought through
> >the matter as carefully as you?[/color]
>
> Considering that the current standard doesn't impose *any* precision on
> sin, I would not expect my DR to be taken seriously.[/color]

So that's your copout for not actually doing anything?
[color=blue]
> Whatever
> considerations led to the current committee decision would be still valid
> after the submission of my DR (no precision specification is *more*
> permissive than *any* precision specification).[/color]

Not even if you were brave enough to pioneer requirements on
precision in the math library? That's not too brave, given that
other languages have done so.
[color=blue]
> If I were to submit a DR, it would be one requiring each implementation
> to document the precision of the 4 arithmetic operations on each floating
> point type. But I know, from past c.s.c discussion, that the committee
> would turn a deaf ear to my DR.[/color]

So that's your copout for not actually doing anything?
[color=blue][color=green]
> >Brave of you. We can call it the Pop Criterion.[/color]
>
> Call it the common sense criterion.[/color]

I won't call it squat if you never submit it as a DR. It's just
more idle carping.

P.J. Plauger
Dinkumware, Ltd.

http://www.dinkumware.com

[Dik T. Winter]

Re: Sine code for ANSI C

In article <c8l06g$34v$18@ sunnews.cern.ch > Dan.Pop@cern.ch (Dan Pop) writes:[color=blue]
> In <HxzJFy.HxG@cwi .nl> "Dik T. Winter" <Dik.Winter@cwi .nl> writes:
>[color=green]
> >To compute the precision of the result you have to do much more than is
> >feasable in most cases.[/color]
>
> What is the value of a result with unknown precision? To repeat my
> earlier example, if you don't know if it's 1 +/- 0.001 or 1 +/- 1000,
> the value of 1 is of pretty little use.[/color]

I have precisely stated what the result means in numerical algebra.
But if you do not trust it see that every use of linear systems with
orders about 40 is abandoned. Or look up in books by Wikinson and
Reinsch. Getting guaranteed results is very expensive. The Karlsuhe
group has even had IBM design specific extensions to their 370 to
make it even feasable.
--
dik t. winter, cwi, kruislaan 413, 1098 sj amsterdam, nederland, +31205924131
home: bovenover 215, 1025 jn amsterdam, nederland; http://www.cwi.nl/~dik/

[Dr Chaos]

Re: Sine code for ANSI C

On Fri, 21 May 2004 14:43:25 GMT, P.J. Plauger <pjp@dinkumware .com> wrote:[color=blue][color=green]
>>
>> For the garbage value, see above. Set errno to any value you consider
>> sensible and document it. ERANGE would be fine with me.[/color]
>
> Okay, so a programmer knows that a call to sine was silly if
>
> 1) the return value is between -1 and +1, and[/color]

I'd think about NaN in an IEEE context.
[color=blue]
> 2) errno indicates that some function value is either too large to
> represent (overflow) or too small to represent (underflow).
>
> That should stand up well to criticism.[/color]

in the case of sin(), yes, because large arguments in floating point
really are mathematically silly. The same is not true of say logarithmic
gamma, or just plain logarithm.

[P.J. Plauger]

Re: Sine code for ANSI C

"Dr Chaos" <mbkennelSPAMBE GONE@NOSPAMyaho o.com> wrote in message
news:slrncb2hbn .ilf.mbkennelSP AMBEGONE@lyapun ov.ucsd.edu...
[color=blue]
> On Fri, 21 May 2004 14:43:25 GMT, P.J. Plauger <pjp@dinkumware .com> wrote:[color=green][color=darkred]
> >>
> >> For the garbage value, see above. Set errno to any value you consider
> >> sensible and document it. ERANGE would be fine with me.[/color]
> >
> > Okay, so a programmer knows that a call to sine was silly if
> >
> > 1) the return value is between -1 and +1, and[/color]
>
> I'd think about NaN in an IEEE context.[/color]

And sin(x) suddenly becomes NaN beyond some arbitrary number of
radians? I'm still waiting for the justification for any particular
magnitude of the cutoff.
[color=blue][color=green]
> > 2) errno indicates that some function value is either too large to
> > represent (overflow) or too small to represent (underflow).
> >
> > That should stand up well to criticism.[/color]
>
> in the case of sin(), yes, because large arguments in floating point
> really are mathematically silly.[/color]

They really aren't, but this seems to be a point too subtle for
many to grasp.
[color=blue]
> The same is not true of say logarithmic
> gamma, or just plain logarithm.[/color]

But the repeated argument is that the silliness stems from the fact
that the floating-point value x stands for the interval (x - 1ulp,
x + 1ulp). If it didn't, it would stand for the value it seems to
stand for, and the library writer might be obliged to actually
compute the function value corresponding to it. Well, what's sauce
for the goose is sauce for the gander. Why should sine get progressively
fuzzier as that interval covers a wider range of function values and
not require/allow/expect exactly the same from every other function
that does the same? The exponential suffers from *exactly* the
same fuzziness problem, to take just the best established of the
fast-moving functions I mentioned earlier. Aren't we doing all our
customers a disservice by giving them a misleading value for exp(x)
when x is large?

P.J. Plauger
Dinkumware, Ltd.

http://www.dinkumware.com

[John Cochran]

Re: Sine code for ANSI C

In article <4tksc.21927$ZQ .18437@nwrddc03 .gnilink.net>,
P.J. Plauger <pjp@dinkumware .com> wrote:
SNIP...
[color=blue]
>
>But the repeated argument is that the silliness stems from the fact
>that the floating-point value x stands for the interval (x - 1ulp,
>x + 1ulp). If it didn't, it would stand for the value it seems to
>stand for, and the library writer might be obliged to actually
>compute the function value corresponding to it. Well, what's sauce
>for the goose is sauce for the gander. Why should sine get progressively
>fuzzier as that interval covers a wider range of function values and
>not require/allow/expect exactly the same from every other function
>that does the same? The exponential suffers from *exactly* the
>same fuzziness problem, to take just the best established of the
>fast-moving functions I mentioned earlier. Aren't we doing all our
>customers a disservice by giving them a misleading value for exp(x)
>when x is large?[/color]

Huh?

What large values for exp(x)?

For most of the floating point implementations out there, an argument
of about 709.78 give or take will give a result close to the upper limit
representable. The main difference for sin() vs exp() is that for all of
the values you can pass to exp(), there is an answer and it is possible
to come up with a floating point number closest to the correct answer to
infinite precision. However, with sin() once you get to an input where the
magnitude of the ulp is greater than 2 pi, it becomes impossible to decide
what the correct answer is.

[Dan Pop]

Re: Sine code for ANSI C

In <hAorc.126$No1. 32@nwrddc02.gni link.net> "P.J. Plauger" <pjp@dinkumware .com> writes:
[color=blue]
>"Dan Pop" <Dan.Pop@cern.c h> wrote in message
>news:c8kvub$34 v$17@sunnews.ce rn.ch...
>[color=green][color=darkred]
>> >: > Assume all the bits of the argument are exact and the *only* source
>> >: > of fuziness is the magnitude of the value. By the time the least
>> >: > significant bit of the mantissa means more than 2 * pi, the fuziness
>> >: > is complete (any result in the -1..1 range is correct).
>> >: >
>> >: > Is that concrete enough for you?
>> >:
>> >: Yes, it's concrete enough to show that you still don't get it.
>> >
>> >So what you're committing to is the requirement that the sine
>> >should *not* be computed for any value where 1 ulp weighs more
>> >than 2*pi.[/color]
>>
>> That's too strong. I would still expect
>>
>> sin(DBL_MAX) * sin(DBL_MAX) + cos(DBL_MAX) * cos(DBL_MAX)
>>
>> to be reasonably close to 1.0, if not equal to it. But, other than that,
>> any value in the [-1, 1] range is OK once 1 ulp exceeds 2*pi.[/color]
>
>There you go again, getting un-concrete. You'll have to specify
>what "reasonably " means, or some nit-picker will take apart
>your revision of the sine specification.[/color]

As close to 1.0 as the C standard requires 2.0 - 1.0 to be ;-)

If you're not interested in a serious discussion, please don't waste my
time.
[color=blue][color=green][color=darkred]
>> >And presumably, sine should deliver the best possible
>> >approximatio n to any argument less than this magnitude.[/color]
>>
>> Not even that. If the value in question covers the range [min, max],
>> sine should provide the best possible approximation of sin(x), where x
>> is one value in the range [min, max].[/color]
>
>Okay, then I see no reason why we shouldn't apply this principle
>to *every* math function. That would mean, for example, that
>the range of permissible values can be huge for large arguments
>to exp, lgamma, tgamma, as well as many combinations of arguments
>to pow.[/color]

Precisely. I didn't mean sine to be treated as an exception.
[color=blue][color=green][color=darkred]
>> >(This
>> >is fine with me, since it still requires an argument reduction
>> >using about 2 1/2 words of precision, so it will break the
>> >vast majority of sine functions that have been in use for the
>> >past several decades.[/color]
>>
>> Only the ones that have been *misused*. If I measure the volume of a
>> sphere with a precision of 1 cubical centimeter, does it make any sense
>> to compute *and display* its radius with picometer precision?[/color]
>
>No. But if you specify the behavior of a math function and don't
>provide a way to specify the uncertainty in its arguments, does
>it make any sense to demand that the math functions know how
>much uncertainty is there?[/color]

The amount of uncertainty assumed by the library is dictated by the
nature of the floating point representation being used. In most real
life applications, the amount is higher than that, but this is none of the
library implementor concern.
[color=blue]
>I've been talking about the requirements
>on the authors of math functions and you keep going back to how
>they should be properly used. Two different universes of discourse.[/color]

They stop being different once common sense is applied. The implementor
need not waste any resources providing more precision than a proper usage
of the function requires. The additional "precision" is garbage, anyway:
only a fool would look for *significant* information in those bits.
[color=blue][color=green][color=darkred]
>> >But hey, we're covered, so I'm happy.)
>> >Of course, you still need to specify what to return for the
>> >garbage value, and what error to report (if any). I await your
>> >further guidance.[/color]
>>
>> For the garbage value, see above. Set errno to any value you consider
>> sensible and document it. ERANGE would be fine with me.[/color]
>
>Okay, so a programmer knows that a call to sine was silly if
>
>1) the return value is between -1 and +1, and
>
>2) errno indicates that some function value is either too large to
>represent (overflow) or too small to represent (underflow).[/color]

The *exact* meaning of ERANGE is as documented for each function using
it. But, as I said, use *any* value *you* consider sensible and document
it.

Again, I'm not interested in this kind of idioting nit picking...
[color=blue][color=green][color=darkred]
>> >So you're willing to stand behind this criterion even after I
>> >described, in a later posting, the various ways in which it was
>> >arbitrary?[/color]
>>
>> Yup. I didn't buy any of your sophistry.
>>[color=darkred]
>> >If this is submitted as a DR to WG14, you're
>> >willing, as it's original author and proposer, to defend it
>> >against cranky attacks by people who haven't thought through
>> >the matter as carefully as you?[/color]
>>
>> Considering that the current standard doesn't impose *any* precision on
>> sin, I would not expect my DR to be taken seriously.[/color]
>
>So that's your copout for not actually doing anything?[/color]

Why should I do *anything*? Are we discussing here the C standard or your
implementation of the C library? Are you so stupid as to be unable to
tell the difference between the two?

I have never expressed any objection to what the C standard says about the
<math.h> stuff, so why the hell should I have to submit *any* DR on this
issue?

Dan
--
Dan Pop
DESY Zeuthen, RZ group
Email: Dan.Pop@ifh.de

[Dr Chaos]

Re: Sine code for ANSI C

On Mon, 24 May 2004 10:51:44 GMT, P.J. Plauger <pjp@dinkumware .com> wrote:[color=blue]
> "Dr Chaos" <mbkennelSPAMBE GONE@NOSPAMyaho o.com> wrote in message
> news:slrncb2hbn .ilf.mbkennelSP AMBEGONE@lyapun ov.ucsd.edu...
>[color=green]
>> On Fri, 21 May 2004 14:43:25 GMT, P.J. Plauger <pjp@dinkumware .com> wrote:[color=darkred]
>> >>
>> >> For the garbage value, see above. Set errno to any value you consider
>> >> sensible and document it. ERANGE would be fine with me.
>> >
>> > Okay, so a programmer knows that a call to sine was silly if
>> >
>> > 1) the return value is between -1 and +1, and[/color]
>>
>> I'd think about NaN in an IEEE context.[/color]
>
> And sin(x) suddenly becomes NaN beyond some arbitrary number of
> radians?[/color]

Yes.
[color=blue]
> I'm still waiting for the justification for any particular
> magnitude of the cutoff.[/color]

When the least significant bit of the argument induces a change of at
least 2 pi its value (so that the output is essentially indeterminate
with a 1 lsb fluctuation), then it is extremely likely that the
programmer is making a blunder and ought to know.
[color=blue]
>[color=green][color=darkred]
>> > 2) errno indicates that some function value is either too large to
>> > represent (overflow) or too small to represent (underflow).
>> >
>> > That should stand up well to criticism.[/color]
>>
>> in the case of sin(), yes, because large arguments in floating point
>> really are mathematically silly.[/color]
>
> They really aren't, but this seems to be a point too subtle for
> many to grasp.[/color]

Explain a situation, using computations with standard fixed-precision
floating point, when computing sin of a value which was passed in
standard floating poitn formats would be valuable for magnitudes
beyond the above.
[color=blue][color=green]
>> The same is not true of say logarithmic
>> gamma, or just plain logarithm.[/color]
>
> But the repeated argument is that the silliness stems from the fact
> that the floating-point value x stands for the interval (x - 1ulp,
> x + 1ulp). If it didn't, it would stand for the value it seems to
> stand for, and the library writer might be obliged to actually
> compute the function value corresponding to it. Well, what's sauce
> for the goose is sauce for the gander. Why should sine get progressively
> fuzzier as that interval covers a wider range of function values and
> not require/allow/expect exactly the same from every other function
> that does the same? The exponential suffers from *exactly* the
> same fuzziness problem,[/color]

It only does in the legalistic sense, but not in the useful esense
because in those cases people would be comparing magnitudes, and
possibly dividing or taking logarithms of the large results.
(And they should have been taking logarithms first most likely)

And the error relative to the value remains small, and that's what
people almost always care about using exp.

At some point you have to return +Inf for exp() of too large an
argument? Is it *really* positive infinity? No. It is far from
positive infinity, well, infinitely so.

There's an arbitrary cutoff right there. Why aren't you forcing
people to convert to an arbitrary precision format?
[color=blue]
> to take just the best established of the
> fast-moving functions I mentioned earlier. Aren't we doing all our
> customers a disservice by giving them a misleading value for exp(x)
> when x is large?[/color]

sin and cos are different from exp, and especially lgamma.

I would be pleased with a library which signaled excessively large arguments
for sin/cos. It's a blunder, like dereferencing NULL.

For a library which had arguments of arbitrary precision, then yes,
it would be a fine idea to get sin() well.
[color=blue]
> P.J. Plauger
> Dinkumware, Ltd.
> http://www.dinkumware.com[/color]