CSR Sparse Matrix Query

Mads Fredrik Skoge Hoel mfhoel at ifi.uio.no
Sun Apr 22 19:36:19 CDT 2007


I was just dealing with this a few moments ago. Remember to make sure the
row and col array start at 0 (and stop at nnz and n-1 respectively). I
just did that mistake.

Regards,
Mads

>
>   You can also use MatCreateSeqAIJWithArrays() and in parallel
> MatCreateMPIAIJWithArrays(). These are intended for those who already
> form the matrix in CSR format.
>
>    Barry
>
>
> On Sat, 21 Apr 2007, Randall Mackie wrote:
>
>> 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.
>> >
>>
>>
>
>





More information about the petsc-users mailing list