Reuse matrix and vector

Jed Brown jed at
Thu Nov 5 03:56:27 CST 2009

jarunan at wrote:
> Hello,
> I would like to reuse matrix and vector to save computing time. As the
> result in last iterations are not similar to the one from my old solver,
> so I am not sure if I program with PETSc the right way, especially,
> resetting values to matrix with MatMPIAIJSetPreallocationCSR().

I suspect you are using a stale preconditioner, but
MatMPIAIJSetPreallocationCSR should not be called every iteration.

> Please take a look, here is how I program it:
> At the beginning of the program I create the vector and matrix.
> call
> VecCreateMPI(PETSC_COMM_WORLD,istorf_no_ovcell,PETSC_DETERMINE,rhs,ierr)
> call VecDuplicate(rhs,sol,ierr)
> call MatCreate(PETSC_COMM_WORLD,Ap,ierr)
> call
> MatSetSizes(Ap,istorf_no_ovcell,istorf_no_ovcell,PETSC_DETERMINE,PETSC_DETERMINE,ierr)
> call MatSetType(Ap,MATMPIAIJ,ierr)
> Then, in each loop I reset values in the vector and the matrix.
> do niter = 1,maxiter
>     call
> VecSetValues(rhs,w,gindex_issu(1:w),f_issu(1:w),INSERT_VALUES,ierr)
>     call VecAssemblyBegin(rhs,ierr)
>     call VecAssemblyEnd(rhs,ierr)
>     call MatMPIAIJSetPreallocationCSR(Ap,rowind,columnind,A,ierr)

It's better to call this once at the beginning and update with
MatSetValues (insert one row at a time, it's very little code).
MatMPIAIJSetPreallocationCSR reallocates memory because it doesn't check
to see if the nonzero pattern has changed (because that's not what it's

>     call MatAssemblyBegin(Ap,MAT_FINAL_ASSEMBLY,ierr)
>     call MatAssemblyEnd(Ap,MAT_FINAL_ASSEMBLY,ierr)
>     call solve_system

Where are you calling KSPSetOperators?  You have to call this every time
the matrix changes.


-------------- next part --------------
A non-text attachment was scrubbed...
Name: signature.asc
Type: application/pgp-signature
Size: 261 bytes
Desc: OpenPGP digital signature
URL: <>

More information about the petsc-users mailing list