Re: Locate an underflow

[Paul van Delst <paul.vandelst(at)noaa.gov>]

William Thompson wrote:
> 
> Paul van Delst <paul.vandelst@noaa.gov> writes:
> 
> >... If, on running said code, I get a crapload of underflow errors, it's an
> >indication that that either a) the code hasn't been tested very well or b) the
> >programmer didn't really think about the problem enough ...
> 
> I disagree.  It's exceedingly easy to get underflow errors, and extremely
> difficult to program around them.  For example, a simple Gaussian
> 
>         Y = A*exp(-((X-X0)/Sig)^2)
> 
> is almost guaranteed to generate underflow errors.  At some point this is going
> to be indistinguishable from zero.  You'd have to jump through hoops to avoid
> getting the completely useless underfloat messages.

Not really. what about something like (assuming double precision):

  tolerance = (MACHAR(/DOUBLE)).EPS

  y = DBLARR( N_ELEMENTS( X ) )
  xarg = ((X-X0)/Sig)^2
  index = WHERE( xarg < tolerance, count )
  IF ( count GT 0 ) THEN $
    y = A*exp(-xarg[index])

Or, as I mentioned in my reply to Craig:

  y = A * gaussian_function( (X-X0)/Sig )

(or similar) which contains all the bits and pieces for checking. Craig also provided a
method of avoiding the underflows.

paulv

-- 
Paul van Delst           A little learning is a dangerous thing;
CIMSS @ NOAA/NCEP        Drink deep, or taste not the Pierian spring;
Ph: (301)763-8000 x7274  There shallow draughts intoxicate the brain,
Fax:(301)763-8545        And drinking largely sobers us again.
                                         Alexander Pope.