How to efficiently change just the diagonal vector in a matrix at every time step

[Sean Dettrick]

One way to do it is to have two Mats, A and B, and a Vec, D, to store  
the diagonal.  A is constructed only on the first step.  On subsequent  
steps, A is copied into B, and then D is added to the diagonal:

   ierr = MatCopy( A, B, SAME_NON_ZERO_PATTERN );
   ierr = MatDiagonalSet( B, D, ADD_VALUES );

The KSP uses B as the matrix, not A.

I don't know if this approach is efficient or not.  Can anybody comment?

Thanks,
Sean



On May 9, 2008, at 2:33 PM, Ben Tay wrote:

> Hi,
> I have a matrix and I inserted all the relevant values during the  
> 1st step. I'll then solve it. For the subsequent steps, I only need  
> to change the diagonal vector of the matrix before solving. I wonder  
> how I can do it efficiently. Of course, the RHS vector also change  
> but I've not included them here.
>
> I set these at the 1st step:
>
> call  
> KSPSetOperators 
> (ksp_semi_x,A_semi_x,A_semi_x,SAME_NONZERO_PATTERN,ierr)
>
> call KSPGetPC(ksp_semi_x,pc_semi_x,ierr)
>
>   ksptype=KSPRICHARDSON
>
>   call KSPSetType(ksp_semi_x,ksptype,ierr)
>
>   ptype = PCILU
>
>   call PCSetType(pc_semi_x,ptype,ierr)
>
>   call KSPSetFromOptions(ksp_semi_x,ierr)
>
>   call KSPSetInitialGuessNonzero(ksp_semi_x,PETSC_TRUE,ierr)
>
>   tol=1.e-5
>
>   call  
> KSPSetTolerances 
> (ksp_semi_x 
> ,tol 
> ,PETSC_DEFAULT_DOUBLE_PRECISION 
> ,PETSC_DEFAULT_DOUBLE_PRECISION,PETSC_DEFAULT_INTEGER,ierr)
>
> and what I did at the subsequent steps is:
>
> do II=1,total
>   call MatSetValues(A_semi_x,1,II,1,II,new_value,INSERT_VALUES,ierr)
>
> end do
>
> call MatAssemblyBegin(A_semi_x,MAT_FINAL_ASSEMBLY,ierr)
>
> call MatAssemblyEnd(A_semi_x,MAT_FINAL_ASSEMBLY,ierr)
>
> call KSPSolve(ksp_semi_x,b_rhs_semi_x,xx_semi_x,ierr)
>
> I realise that the answers are slightly different as compared to  
> calling all the options such as KSPSetType, KSPSetFromOptions,  
> KSPSetTolerances at every time step. Should that be so? Is this the  
> best way?
>
> Also, I can let the matrix be equal at every time step by fixing the  
> delta_time. However, it may give stability problems. I wonder how  
> expensive is these type of value changing and assembly for a matrix?
>
> Thank you very much.
>
> Regards.
>