Actual source code: ex262f.F90

  1:       program main
  2: #include <petsc/finclude/petscmat.h>
  3:       use petscmat
  4:       implicit none

  6:       Mat      A,B
  7:       PetscErrorCode ierr
  8:       PetscScalar, pointer :: km(:,:)
  9:       PetscInt three,one
 10:       PetscInt idxm(1),idxmj(1),i,j
 11:       PetscMPIInt rank,size

 13:       PetscCallA(PetscInitialize(ierr))

 15:       PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr))
 16:       PetscCallMPIA(MPI_Comm_size(PETSC_COMM_WORLD,size,ierr))

 18:       PetscCallA(MatCreate(PETSC_COMM_WORLD,A,ierr))
 19:       three = 3
 20:       PetscCallA(MatSetSizes(A,three,three,PETSC_DECIDE,PETSC_DECIDE,ierr))
 21:       PetscCallA(MatSetBlockSize(A,three,ierr))
 22:       PetscCallA(MatSetUp(A,ierr))
 23:       PetscCallA( MatDuplicate(A,MAT_DO_NOT_COPY_VALUES,B,ierr))
 24:       one = 1
 25:       idxm(1) = 0
 26:       allocate (km(three,three))
 27:       do i=1,3
 28:         do j=1,3
 29:            km(1,1) = i + j
 30:            idxm(1) = i - 1 + 3*rank
 31:            idxmj(1) = j - 1 + 3*rank
 32:            PetscCallA(MatSetValues(B, one, idxm, one, idxmj, reshape(km, [three*three]), ADD_VALUES, ierr))
 33:         enddo
 34:       enddo

 36:       PetscCallA(MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY,ierr))
 37:       PetscCallA(MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY,ierr))
 38:       PetscCallA(MatView(B,PETSC_VIEWER_STDOUT_WORLD,ierr))

 40:       PetscCallA(MatDestroy(A,ierr))
 41:       PetscCallA(MatDestroy(B,ierr))

 43:       deallocate(km)
 44:       PetscCallA(PetscFinalize(ierr))
 45:       end

 47: !/*TEST
 48: !
 49: !   test:
 50: !     nsize: 2
 51: !
 52: !TEST*/