Passing function handles and inline functions to fsolve, and family

[David Bateman]

Dear All,

In the endless goal of further compatiability I've been trying to get all of
the functions of octave that take a function as an argument (generally as
a string) to also allow function handles and inline functions to be accepted.
These functions are

daspk
dasrt
dassl
fsolve
lsode
odessa
quad

The function "quad" is not problematic, however the others are.
There are two basic problems, firstly the multiple return args of the
inline functions and secondly the differing interfaces between these
functions

Compatiability with the matlab fsolve function means that the inline
or function handle that is passed should be a function that returns
multiple arguments for the function evaluation, jacobian, etc. This is
problematic since the existing code treats these as seperate functions.
In fact the underlying fortran code also treats these as seperate functions.

The brain dead approach to allowing the inline function to return
multiple arguments is the wrap the function so that only the first,
second or third return argument is used in the various functions. This
change is relatively simple but comes at the cost of upto three times
as many function evaluations than are really needed.

As the underlying fortran code requires seperate functions, and
barring altering this fortran code or replacing it, we are left with
the only other solution of caching the return function evaluations,
jacobians, etc and using them on the next call. This requires more
extensive modifications but is marginally workable. The presence of
the jacobian as a second return value is flagged in fsolve_ops in a
similar manner to matlab.

The major problem is the DASRT and DASSL functions where the call for
function evaluation from the underlying fortran code doesn't contain
all of the required arguments to allow the jacobian to be
calculated. Barring altering the underlying fortran code my only
solution is to pass an illegal value for the missing value in such a
way as to invalidate the use of the cached jacobian values. A similar
problem exists for the ODESSA inhomogeneity matrix.

Another minor problem with this approach is that the cached values are
at the moment stored as static variables. So the code as it stands is
not re-entrant. This means that if we can use gfortran to compile the
stuff in libcruft and allow recursive calls, then we'll still have
problems with recursion due to this non re-entrant code. Writing as 
re-entrant is not that difficult and so I haven't treated this issue
yet.

My other problem with the consistence of the interface to the
functions is with odessa and dasrt with also allow a constraint
function or inhomogeneity matrix to be given. In this case the dasrt
function currently has the function evaluation and jacobian passed as
one argument and the constraint function passed as another, while
odessa has all three functions passed as a single argument.  My gut
feeling is that the odessa interface is the right way to go, but I
hesitate to change the dasrt interface due to backward compatiability
issues. Either I'd have to drop backward compatiability or keep a
large amount of cruft for the purpose of dealing with it.

In any case I've written a complete patch to deal with this and it is
attached here. It seems to work fine on the test cases I've given it.
But due to the issues above I'm not entirely happy with it. I'd be happy
if users of these functions would suggest better alternatives to this
patch, but barring that I'd suggest this patch be applied...

Regards
David

-- 
David Bateman                                address@hidden
Motorola CRM                                 +33 1 69 35 48 04 (Ph) 
Parc Les Algorithmes, Commune de St Aubin    +33 1 69 35 77 01 (Fax) 
91193 Gif-Sur-Yvette FRANCE

The information contained in this communication has been classified as: 

[x] General Business Information 
[ ] Motorola Internal Use Only 
[ ] Motorola Confidential Proprietary

-- 
David Bateman                                address@hidden
Motorola CRM                                 +33 1 69 35 48 04 (Ph) 
Parc Les Algorithmes, Commune de St Aubin    +33 1 69 35 77 01 (Fax) 
91193 Gif-Sur-Yvette FRANCE

The information contained in this communication has been classified as: 

[x] General Business Information 
[ ] Motorola Internal Use Only 
[ ] Motorola Confidential Proprietary

changelog-handle
Description: Text document

patch.handle
Description: Text document