Re: matrix type caching, triangualar and Cholesky solvers

[David Bateman]

John W. Eaton wrote:

>On 27-Apr-2006, David Bateman wrote:
>
>| Here is a patch that allows whether a full matrix is upper, lower or
>| positive definite to be probed, cached and the method of solving linear
>| equations with this matrix with the "/" or "\" operators automatically
>| chosen based on this information. This seems to be about 40% for
>| positive definite matrices and about 100 times for triangular matrices
>| using atlas on my machine. I tried to test all possible cases of the use
>| of this code before submitting it, but as there is no testing code
>| whatsoever for the "\" and "/" full operators yet, it wasn't easy to
>| propose appropriately adapted test code in the absence existing test code.
>| 
>| This should supersede the chol.oct file from octave-forge and the
>| ov-tri-mat type as well.. John do you want it committed?
>
>I'd say yes, but I have a couple of questions.
>
>In ltsolve, you have
>
>  if (typ == MatrixType::Permuted_Lower)
>    {
>      (*current_liboctave_error_handler)
>         ("Permuted triangular matrix not implemented");
>    }
>  
>
Matlab allows permuted triangular matrices to be treated, and we
probably should as well. Then if you do something like "[l,u] = lu(a); x
= l\u\b;" The permuted lower triangular matrix will be accelerated as
well. I think matlab uses the dmperm code on dense matrices to decide if
the matrix is a permuted triangular matrix, but I'm not sure. In any
case I have some of the code to allow permuted triangular matrices
already in place for a later fix.

Note matlab also treated full Hessenberg seperately as well and I have
treated this case either yet.

>Does the current code handle this case or is this something that would
>fail now?  I guess there is no way to actually generate this case
>(other than manually tagging a matrix), so it is OK?
>  
>
[l,u] = lu(a) definite generates a permuted lower triangular matrix..

>What is the intended usage of the MatrixType::info function?  Why does
>it use the warning handler?
>  
>
This is a reflection of the same thing for sparse matrices, where the
spparams("spumoni",1) flag allows printing of what the code is doing
during the solve. This is not fully active yet as I'm not sure what flag
to use for the print of matrix solve information.

One other point is that the MatrixType and SparseType should probably be
combined, and the sparse_type and matrix_type methods in ov-base-mat.h
and ov-base-sparse.h moved to ov.h as a single method. This will
simplify code in matric_type.cc and avoid the dynamic_cast..

>| Does someone else want to try it out?
>
>That might be good too.
>
>jwe
>
>  
>
Cheers
David

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

The information contained in this communication has been classified as: 

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