Re: [Lapackpp-devel] LaBandFactDouble Question

[Christian S. <sti...@tu...>]

Dear Christophe,

Christophe Picard schrieb:
> I am trying to use the banded matrices to solve pde, and I notice that the
> return arguments in the function LaBandMatFactorize are weird. 

Absolutely. Yes, in fact *all* the *FactDouble classes are pretty much 
messed up. In gfd.h, I added a comment which made this clear, but in the 
other Fact headers I forgot to add such a comment (will be included 
immediately).

In other words, the design of these classes is screwed anyway. If you 
found a way to make them useable, then you would be very welcome to 
contribute your code to lapackpp. I will happily throw out the wrong 
code and replace it by your code.

> I mean that
> the pivot and upper/lower matrix should be save in the member B and not in
> the original matrix:

Err, yes, but the whole class needs more thought. It depends on how you 
construct the LaBandFactDouble object. If you use the constructor 
LaBandFactDouble(LaBanMatDouble& A), then the member variable B 
references the same memory as A and the function LaBandMatFactorize will 
use the matrix A in-place, regardless whether A or B is passed as parameter.

I'm not sure how such a factorization is used in practical numerics. It 
might well be possible that you should throw away the whole 
implementation of bfd.h and invent your own. On the other hand, I 
recently added the class LaGenQRFactComplex in gfqrc.h and I think its 
design is one possibility to efficiently model matrix factorization as a 
class. Maybe you can try to copy the implementation model of gfqrc.h 
into one that uses LaBandMatDouble, but that's just one possibility out 
of many.

> So my question is in bfd.h : (...)
> should it be replace with something more along this (even maybe overload
> the function)
> 
> 00138 inline void LaBandMatFactorize(LaBandMatDouble &A, LaBandFactDouble
> &AF)
> 00139 {
> 00140     integer n = A.size(1), m = n, LDA = A.gdim(0);
> 00141     integer KL = A.subdiags(), KU = A.superdiags(), info=0;
> 00142
> 00143     F77NAME(dgbtrf)(&m, &n, &KL, &KU, &AF.B(0,0), &LDA,
> &(AF.pivot()(0)), &info);
> 
> otherwise the factorization is not saved into AF but in A, making AF.B
> useless, and the function LaLinearSolve unusable.

Absolutely. Throw out the old errornous function and replace it by your 
own function, and I would happily include that in the official lapackpp 
repository. Especially if you have tested and verified your implementation.

> This is just a question/remark.

Thank you very much for this feedback.

Regards,

Christian Stimming