CSR Sparse Matrix Query

Randall Mackie randy at geosystem.us
Sat Apr 21 08:53:28 CDT 2007


Hi Tim,

There's no need to generate a dense matrix if you're dealing with sparse
matrices. There are many approaches, but if you just want to get your
existing matrix into PETSc, you can use something like the following (written
in Fortran, but almost the same in C):

Assume the global vector length is np:

       call VecCreateMPI(PETSC_COMM_WORLD,PETSC_DECIDE,np,b,ierr)
       call VecDuplicate(b,xsol,ierr)
       call VecGetLocalSize(b,mloc,ierr)
       call VecGetOwnershipRange(b,Istart,Iend,ierr)

       do i=Istart+1,Iend
         loc(i)=i-1
       end do


       ALLOCATE (d_nnz(mloc), o_nnz(mloc), STAT=istat)

Note: the d_nnz and o_nnz are used to determine the pre-allocation of the global matrix
according to the PETSc manual (basically just looping through and figuring out which
elements are in d_nnz and which are in o_nnz).

Then you can create the matrix:


       call MatCreateMPIAIJ(PETSC_COMM_WORLD,mloc,mloc,np,np,
      .     PETSC_NULL_INTEGER, d_nnz, PETSC_NULL_INTEGER,
      .     o_nnz,P,ierr)


And then I had a subroutine to fill it like this (the critical part is the check
whether that part of the matrix is owned by the current processor):

       jj=0

       do i=1,l
         do k=2,n
           do j=2,m

             jj=jj+1
             row = jj-1

             IF (jj >= Istart+1 .and. jj <= Iend) THEN

               if(k.ne.2) then
                 fact=-rijcy*dx*dy/z(k-1)       ! Hxijc
                 ic=ic+1
                 v(ic)=fact
                 col(ic)=ijkhx(i,j,k-1)-1
               end if

[snip code]

               call MatSetValues(A,i1,row,ic,col,v,INSERT_VALUES,
      .             ierr)


[and finally at the end you assemble the matrix]

       call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)
       call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)


Since these early days, I've switched over to using Distributed Arrays,
which in my opinion are much better and much easier to use.

Hope this helps,

Randy


Dr. Timothy Stitt wrote:
> Hi all,
> 
> I have just begun to use PETSc and I have already ran into a problem. 
> Apologies if this has been dealt with before.
> 
> My code currently generates its own non-distributed Hamiltonian sparse 
> matrix in CSR format (due to the large size of the dense 
> representation). I want to exploit PETSc's (actually SLEPc's) 
> eigensolver routines and hence I need to convert my CSR representation 
> into PETSc CSR representation for use by the eigensolver.
>  From reading the documentation I am under the impression that you 
> initially supply the dense matrix to the PETSc matrix routines for 
> conversion to PETSc CSR representation. In my case I do not generate the 
> dense matrix. I would be grateful if someone could guide me on to how to 
> create a PETSC sparse matrix from my own in-code generated sparse matrix.
> 
> Thanks in advance for any assistance that can be provided.
> 
> Regards,
> 
> Tim.
> 

-- 
Randall Mackie
GSY-USA, Inc.
PMB# 643
2261 Market St.,
San Francisco, CA 94114-1600
Tel (415) 469-8649
Fax (415) 469-5044

California Registered Geophysicist
License No. GP 1034




More information about the petsc-users mailing list