(2nd time) tweaks to erff() threshold values

[Steve Kargl]

On Wed, Aug 28, 2013 at 09:31:21PM +1000, Bruce Evans wrote:
> On Tue, 27 Aug 2013, Steve Kargl wrote:
> 
> > On Tue, Aug 27, 2013 at 09:24:39PM +1000, Bruce Evans wrote:
> >> On Mon, 26 Aug 2013, Steve Kargl wrote:
> 
> >>> ...
> >>> /*
> >>> - * Coefficients for approximation to  erfc in [1.25,1/0.35]
> >>> + *  Domain [1.25,1/0.35], range [-7.043e-10,7.457e-10]:
> >>> + *  |log(x*erfc(x)) + x**2 + 0.5625 - r(x)/s(x)| < 2**-30
> >>>  */
> >>
> >> It is even less clear than for the others that the scaling of the error
> >> is correct.  The comment could make this clearer by giving the ratio of
> >> the correct error scale with the one used in the approximation, at the
> >> endpoints (max and min for this, usually at the endpoints).
> >
> > I'm not sure what you mean here.  In looking at error curves,
> > it turns out that one of the 2 limits in the range is taken
> > from an endpoint.
> 
> I mean that these are absolute errors for a function related to the
> result that is especially obscure in this case -- the function is
> log(x*result).  The error that needs to be minimized is the relative
> error of the result, where "relative" is not quite the ratio of the
> error and the result -- it is the ratio of the error and [result
> rounded to a nearby power of 2; I'm not sure if the rounding is
> up or down].  For the absolute error to be close enough to the
> relative error, the relative error must not vary much across the
> interval, AND the scaling for the abslute error must be compatible.
> This is far from clear here.  I think taking logs changes the scale
> for the absolute error to logarithmic, so this error is not really
> absolute.  Hopefully the scale for the result is similarly logarithmic,
> so that the scales are compatible.  There is also a factor of x before
> taking logs, so the scaling is not logarithmic either.

Before I go down the wrong road, I thought I would check if I'm 
interpreting the above correctly.  To be concrete, in the interval
[1/0.35, 11], we consider the function

(1)   f(x) = log(x * erfc(x)) + x**2 + 0.5625 

and I've minimized the maximum error defined by

(2)   e(x) = f(x) - p(x) / (1 + q(x))

where p(x) = p0 + p1 * z + ..., q(x) = q1 * z + q2 * z**2 + ..., and
z = 1 / x**2.

In my current code, I get the following error curve:

 4e-12 ++-----+-------+------+-------+------+-------+------+-------+-----++
        +      +       +      +       +   ****       +      +       +      *
  3e-12 ++             **                **  **                           **
        |     * **     ***              **    **                          *|
        |     * **    ** *             **      **                        **|
  2e-12 ++    * **    *  **            *        **                       *++
        |     * **    *   *           **         **                      * |
  1e-12 ++    * **    *   *           *           **                    **++
        |     * **   **   **         **            *                    *  |
        |     *****  *     *         *             **                  **  |
      0 ++     *  *  *     *        **              **                **  ++
        |      *  *  *     **       *                *                *    |
 -1e-12 ++     *  *  *      *      **                **              **   ++
        |      *  * **      *      *                  **             *     |
 -2e-12 ++     *  * *       **    **                   **           **    ++
        |      *  ***        *   **                     **         **      |
        |      *   **        **  *                       **       **       |
 -3e-12 ++                    ****                        **    ***       ++
        +      +       +      +       +      +       +     ******   +      +
 -4e-12 ++-----+-------+------+-------+------+-------+------+-------+-----++
        2      3       4      5       6      7       8      9       10     11


If I rewrite Eq. (1) to isolate erfc(x), I arrive at

(3) erfc(x) = F(x) = exp(p(x) / (1+q(x)) - x**2 - 0.5625) / x

and if I understand your comment above, I should define a new
error curve as

(4) E(x) = (erfc(x) - F(x)) / erfc(x)

where erfc(x) is the assumed exact result and E(x) is a relative error.
The relative error curve, using coefficients determined from (2), then
looks like

  4e-12 ++-----+-------+------+-------+------+-------+------+-------+-----++
        +      +       +      +       +   ****       +      +       +      *
  3e-12 ++             **                **  **                           **
        |     * **     ***              **    **                          *|
        |     * **    ** *             **      **                        **|
  2e-12 ++    * **    *  **            *        **                       *++
        |     * **    *   *           **         **                      * |
  1e-12 ++    * **    *   *           *           **                    **++
        |     * **   **   **         **            *                    *  |
        |     *****  *     *         *             **                  **  |
      0 ++     *  *  *     *        **              **                **  ++
        |      *  *  *     **       *                *                *    |
 -1e-12 ++     *  *  *      *      **                **              **   ++
        |      *  * **      *      *                  **             *     |
 -2e-12 ++     *  * *       **    **                   **           **    ++
        |      *  ***        *   **                     **         **      |
        |      *   **        **  *                       **       **       |
 -3e-12 ++                    ****                        **    ***       ++
        +      +       +      +       +      +       +     ******   +      +
 -4e-12 ++-----+-------+------+-------+------+-------+------+-------+-----++
        2      3       4      5       6      7       8      9       10     11

In fact, if I plot e(x) and E(x) on the same (high resolution) figure,
the curves are indistingiushable.  So, is (4) not the right relative
error curve?
 
-- 
Steve