bessel.lisp error codes and noisy quad_qagi



On Nov. 11, 2011, Raymond Toy wrote:
------------
>> quad_qagi.  The rest of the zbesk ierr = stuff is not needed,
>> nor is the final ***MESSAGE FROM ROUTINE DQAGI IN LIBRARY SLATEC.
>> For my purposes, this is noise.
>>
>
>I was thinking of slightly bigger picture than your nint.  Right now,
>that
>message is the only way to know that zbesk might not be doing what you
>think it is.  There ought to be a way for a program to determine that.
>
>For you particular issue, I also see no reason for quad_qagi to continue
>doing any more computations with bogus results from zbesk.  If zbesk
>could
>indicate failure, then quad_qagi could immediately return failure.  That
>would be more robust for your purposes.  It's hard to tell if quad_qagi
>would return success even if zbesk was completely failing.  Worst case:
>zbesk is totally broken and returns 0 for everything.  quad_qagi would
>think everything was fine and give a nice "reliable" answer.
>
----------------------
I agree that a wider and deeper approach would be valuable.

A possible approach to quiet mode error handling
    with Maxima special functions is sketched out
  below:

We need:

1.  a global list,   spfun_errL,   by default empty

2.  a global flag,  spfun_quiet_mode,   by default false

When a layered package program (such as nint)
calls quad_qagi, etc, the calling program first
sets spfun_quiet_mode:true and makes sure that
spfun_errL is an empty list (length zero).

Inside quad_qagi, for example,
between each successive refinement of the value
of the numerical integral (increase in number
of mesh points, etc.) quad_qagi examines the
global list spfun_errL to see if it is empty
 or has length greater than zero, in which
case quad_qagi, etc, returns with an error
code of 6.

(See below for how spfun_errL gets populated
with sublists.)


Although in principle, the action taken by a
program which calls special functions
could vary, depending on the contents
of spfun_errL, quad_qagi could just return a list
with the approximation reached at the previous
iteration, except that an error code = 6 would be
used as the fourth element of the list returned.

If the program calling quad_qagi is nint, that
error code = 6 return from quad_qagi would
trigger an appropriate message to the user of the
program, and nint would be responsible for
resetting the list spfun_errL to an empty
list before further calls, and also setting
the flag spfun_quiet_mode to its default
value false before exiting the program.

Populating spfun_errL:

Each of the Maxima special functions is given an
id number (bessel-j,1) (bessel-y,2), (bessel-i,3),
 (bessel-k,4),   etc in their code.

When a Maxima special function (such as bessel-k) detects an error code
such as is returned from zbesk.lisp, and if spfun_quiet_mode
is true, then instead of printing an error message to
the console, bessel-k conses the sublist [id-num, 
bessel-fun-order,bessel-fun-arg,ierr]
into spfun_errL and returns something appropriate
(perhaps 0.0)  to quad_qagi.

If the integrand needs calls to more than one
special function, and if more than one such call results
in an error, then more than one sublist will appear in
spfun_errL.

An interactive user of numerical values of Bessel functions
would have the option of depending on the usual printed
error warning, or of using the same spfun_quiet_mode.
In the latter case, the populated list spfun_errL could be
used for diagnostics, including the improvement of
the bessel function routines.

If the list spfun_errL is actually initially a Lisp list,
then the corresponding Maxima list could be created via a
plain words translation of the id-num value and ierr value
into, for example, bessel_k, overflow or underflow.

If all Maxima special functions adopted a common meaning
for numerical error return numbers, the job would be
simplified.

Ted