MatLUFactor for Complex Matrices - Returns inverse of the diagonal in LU??

[Paul Dostert]

Barry,
Thank you so much for getting back to me so quickly. Just out of curiosity,
is there a slightly more in depth version of the documentation that I would
be able to look at if I have some other "simple" questions like this in the
future?


On Wed, May 27, 2009 at 5:59 PM, Barry Smith <bsmith at mcs.anl.gov> wrote:

>
>   Paul,
>
>     In PETSc for the AIJ sparse matrix storage format we store the inverse
> of the diagonal entries of U. This is because then when it comes to the
> triangular solves they only involve floating point multiplies and adds (no
> time consuming divisions). I think this is pretty common for sparse matrix
> solvers. You could argue that the MatView() should invert those values so
> what gets printed out is clearer; I view the MatView() for factored matrices
> as only a useful tool for debugging for developers so it just prints out the
> stored values. Yes it is a bit confusing.
>
>    Barry
>
>   For BAIJ matrices it gets even more confusing. With BAIJ the
> factorizations use the little blocks of the matrix as the basic
> computational unit of the factorization. Here we actually compute the exact
> inverse of the diagonal block and store that (a block version of what we do
> for AIJ). Again this is to make the solves as efficient as possible.
>
>
>
>
> On May 27, 2009, at 6:20 PM, Paul Dostert wrote:
>
>  I'm a beginner with PETSc, so please forgive me if this is obvious, but I
>> couldn't seem to find any help in the archives.
>>
>> I'm trying to just the hang of the software, so I've been messing around
>> with routines. I'm going to need complex matrices (Maxwell's with PML) so
>> everything is configured for this. I'm messing around with some very
>> simple
>> test cases, and have a symmetric (but not Hermitian) complex matrix with
>> 2-2i on the diagonal and -1+i on the upper and lower diagonals. I am
>> reading
>> this in from a petsc binary file (again, for testing purposes, eventually
>> I'm going to be just reading in my matrix and RHS). I view the matrix, and
>> it has been read in correctly. I perform LU factorization by doing the
>> following (where ISCreateStride(PETSC_COMM_WORLD,m,0,1,&perm); has been
>> called earlier):
>>
>> MatConvert(A,MATSAME,MAT_INITIAL_MATRIX,&LU);
>> MatFactorInfoInitialize(&luinfo);
>> MatLUFactor(LU,perm,perm,&luinfo);
>>
>> I get that LU is:
>>
>>  (1,1)      0.2500 + 0.2500i
>>  (2,1)     -0.5000
>>  (1,2)     -1.0000 + 1.0000i
>>  (2,2)      0.3333 + 0.3333i
>>  (3,2)     -0.6667
>>  (2,3)     -1.0000 + 1.0000i
>>  (3,3)      0.3750 + 0.3750i
>>  (4,3)     -0.7500
>>  (3,4)     -1.0000 + 1.0000i
>>  (4,4)      0.4000 + 0.4000i
>>  (5,4)     -0.8000
>>  (4,5)     -1.0000 + 1.0000i
>>  (5,5)      0.4167 + 0.4167i
>>  (6,5)     -0.8333
>>  (5,6)     -1.0000 + 1.0000i
>>  (6,6)      0.4286 + 0.4286i
>>
>> Now, I am interpreting this as L being unit on the diagonal and the
>> lower diagonal portion of this "LU" matrix, while U being the diagona
>> + upper of this "LU" matrix. I can interpret this the other way around
>> as well, and it doesn't matter.
>>
>> However, knowing the LU factorization, it is VERY clear the the proper
>> LU decomposition would have the inverse of the diagonal elements
>> presented here. So I believe I should have  LU as:
>>
>>  (1,1)      2.0000 - 2.0000i
>>  (2,1)     -0.5000
>>  (1,2)     -1.0000 + 1.0000i
>>  (2,2)      1.5000 - 1.5000i
>>  (3,2)     -0.6667
>>  (2,3)     -1.0000 + 1.0000i
>>  (3,3)      1.3333 - 1.3333i
>>  (4,3)     -0.7500
>>  (3,4)     -1.0000 + 1.0000i
>>  (4,4)      1.2500 - 1.2500i
>>  (5,4)     -0.8000
>>  (4,5)     -1.0000 + 1.0000i
>>  (5,5)      1.2000 - 1.2000i
>>  (6,5)     -0.8333
>>  (5,6)     -1.0000 + 1.0000i
>>  (6,6)      1.1667 - 1.1667i
>>
>> Is there some reason this returns the inverse of the diagonal entries,
>> or am I completely missing something? Is returning the inverse
>> something standard??
>>
>> Since I'm new, I'm not quite sure where to look for actual source
>> code. Is there a location where the LU factorization code is written
>> and well commented?
>>
>> Thank you very much!
>>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20090527/ae96e33e/attachment.htm>