Re: fsolve improvements

[Jaroslav Hajek]

On Mon, Jan 26, 2009 at 2:31 PM, Carlo de Falco <address@hidden> wrote:
>
> On 26/gen/09, at 10:53, Jaroslav Hajek wrote:
>
>> hello,
>>
>> I have probably finished all the fsolve improvements I intended to do
>> for 3.2. I intend to do some extensive testing in the weeks to come,
>> so changes may yet
>> be made, but the theory is all there.
>>
>> Here's a summary for users interested in solving systems of nonlinear
>> equations:
>>
>> fsolve is now an m-file function, the old minpack code is gone. thus,
>> fsolve can now be called recursively, and will automatically compute
>> in single precision if given single-precision inputs. General
>> rectangular systems should now be supported, solved in least-squares /
>> minimum-norm sense. The interface and options supported via the
>> optimset are mostly matlab-compatible, as well as the return values
>> (where applicable).
>> Essentially the same method as in minpack is used - double dogleg
>> trust region steps with Broyden updating of the jacobian. For large
>> enough (no. of variables > 10) square and overdetermined(!) systems,
>> the QR factorization of the jacobian is updated rather than the
>> jacobian itself, avoiding costly re-factorizations. Broyden updating
>> may be switched off by the user, yielding plain levenberg-marquardt
>> method.
>> To fully benefit from the efficient updating, you should link Octave
>> with the qrupdate library.
>>
>> In particular, fsolve should now be capable of general nonlinear
>> fitting with expected small residuals (see the last test of fsolve).
>> I encourage users to test the new fsolve and report
>> problems/suggestions to this list.
>>
>> cheers
>
> Hi,
> Thanks for the great job!
> I was wondering, does fsolve get any benefit is the jacobian is sparse?
> I'm using 3.0 at the moment but I would have a few nice applications with
> nonlinear PDEs that could be interesting tests for fsolve, so I'd be tempted
> to switch to the development sources if
> sparse jacobians are handled.
> Carlo
>
>

Hi Carlo,

currently, the best strategy if sparse jacobian is supplied by the
user is to disable the Broyden updates, because they destroy the
sparsity and treat the jacobian as dense. In fact, I think I should
make this automatic if sparse jacobian is detected - thanks! As a
result, a sparse system will be solved in each step, otherwise the
code should work normally.
To avoid solving a system in each step, various strategies can be
used, such as Schubert's updating. The limitation here is that the
qrupdate library does not (yet) support sparse QR factorizations. This
can change in the future, depending on my motivation to implement it
(currently, I have little business with sparse systems).

-- 
RNDr. Jaroslav Hajek
computing expert
Aeronautical Research and Test Institute (VZLU)
Prague, Czech Republic
url: www.highegg.matfyz.cz