inverting independent, small dense matrices

[Shao-Ching Huang]

Barry:

How do I determine the "shift" value, to be passed into

   Kernel_A_gets_inverse_A_4(MatScalar *a,PetscReal shift)

?

I know that my 4x4 matrices are non-singular for sure.

Thanks.

Shao-Ching


Barry Smith wrote:
> 
> On Feb 13, 2009, at 3:01 PM, Shao-Ching Huang wrote:
> 
>> Hi
>>
>> I need to invert a bunch of independent, dense 4x4 matrices on each
>> MPI process. (This is separate from the main computation, a parallel
>> Ax=b solve.) I am thinking of using LAPACK's DGESV to invert these 4x4
>> matrices.
>>
> 
>      You would not want to use LAPACK's solvers for such tiny matrices. 
> LAPACK is optimized
> for large dense matrices; the overhead is high for small matrices so 
> that actual performance will
> be very slow calling LAPACK.
> 
>> In my code, is it necessary to call LAPACK directly? Or does PETSc
>> already have a convenient interface to LAPACK?
>>
>      For PETSc's own internal use we have written custom inverters for 
> small matrices (these are used by
> the factorization routines for BAIJ matrices).  They are called 
> Kernel_A_gets_inverse_A_* where * is the block size.
> The prototypes for the functions are in src/inline/ilu.h and the actual 
> source code in files in src/mat/impls/baij/seq/dgef*.c
> 
>    Barry
> 
> 
>> The size 4x4 may be increased to 8x8 or so later, so I am not trying
>> to use the closed-form inversion formula.
>>
>> Thanks,
>>
>> Shao-Ching
>