Aleksandar Donev writes:
> It is simply a
> quadratic equation that needs to be solved, and the answer is C /
> (-B+Sqrt(B**2-A*C)) if the two spheres collide (i.e. if this result is
> nonnegative and exists). I need to be very careful that this gives a
> very accurate answer...
In that case, perhaps a different formula would be appropriate. That
one is notorious for being numerically poor in some cases - sounds
like the kinds of cases you have. I haven't needed highly accurate
quadratic equation solutions myself, so I don't recall the details,
but I'm sure I've seen people discuss more accurate ways to compute
this in the numerically sensitive cases. I could probably even derive
some with a bit of work....but don't think I'll take the time to try
right now.
> However, as NAG is telling me, suspicious IEEE operations are happening:
>
> Warning: Floating divide by zero occurred
> Warning: Floating underflow occurred
....
suspicious indeed. If you are getting things like that, you are down
in the area where, yes, you will be affected by subtle things like
whether intermediate results are kept in higher-precision registers or
flushed to memory. This in turn can be affected by all kinds of
things, such as optimization level, or the presense of diagnostic
output. Nasty.
> How should I write something like this in Fortran 90? Or is using the
> IEEE TR module and its facilities the only way to go?
I'd recommend a less sensitive formula. But as I was typing this,
out shipping guy came in with something that out-prioritizes this,
so I'll beg out on finding or figuring out what the better formulae
might be.
--
Richard Maine | Good judgment comes from experience;
[log in to unmask] | experience comes from bad judgment.
| -- Mark Twain
|