# A^TA + B?

Zhifeng Sheng z.sheng at ewi.tudelft.nl
Fri Oct 26 11:07:01 CDT 2007

```Hong Zhang wrote:

>
> Zhifeng,
>
>>
>> I would like to implement B = A^TA + B.
>>
>> And the transpose multiplication can be done with the following
>> function,
>>
>> PetscErrorCode PETSCMAT_DLLEXPORT MatMatMultTranspose(Mat A,Mat
>> B,MatReuse scall,PetscReal fill,Mat *C) since I would like to save
>> some memory, I don't want to create an intermediate matrix. In the
>> above function, the discription of
>>
>> the argument "scall" is not clear, I know that MAT_INITIAL_MATRIX
>> will create a new empty matrix and the result will be stored.
>
>
> Yes.
>
>>
>> what about MAT_REUSE_MATRIX ? does it mean that the C matrix will be
>
>
> When the product matrix C=A^T*B is previously computed,
> and you want repeating C=A^T*B, in which A and B maintain
> the same non-zero sparse pattern but with different numerical values,
> using
> MAT_REUSE_MATRIX will skip the symbolic computation of C=A^T*B
> and reuse the exiting matrix data structure and memory of C.
>
> Why do you need B = A^TA + B?
> If your matrices are dense, you should use lapack or scalapack.
> Sparse matrix product should be avoided because
> the product matrix C usually is much denser than A and B,
> and C=A^T*B cannot be implemented efficiently in general.
> We provide MatMatMultTranspose() mainly for multigrid computation
> in which the matrices have special data structure.
>
> Best,
>
> Hong
>

I need to implement a least-squares method, in which A^T*A will be
computed, however, since I use finite elements, local matrices will be
computed and then added to the global matrix. The A matrix is very
sparse, e.g 9 non-zeros in a row.

Actually, in this case, I only need to do b = a^T*a + b on local
matrices, which are very small e.g (12 by 12). The elementary matrices
will be assembled in the end. Should I use lapack in combination with Petsc?

Best regards
Zhifeng

```