1program main 2#include <petsc/finclude/petscvec.h> 3#include <petsc/finclude/petscmat.h> 4 5use petscvec 6use petscmat 7 8implicit none 9 10 Mat A 11 PetscInt,parameter :: n=5,m=5 12 PetscScalar,parameter :: two =2.0, one = 1.0 13 PetscInt,pointer,dimension(:) :: dnnz,onnz 14 PetscInt :: i,rstart,rend,M1,N1 15 PetscErrorCode ierr 16 17 PetscCallA(PetscInitialize(ierr)) 18 19 allocate(dnnz(0:m-1)) 20 allocate(onnz(0:m-1)) 21 22 do i=0,m-1 23 dnnz(i) = 1 24 onnz(i) = 1 25 end do 26 27 PetscCallA(MatCreateAIJ(PETSC_COMM_WORLD,m,n,PETSC_DETERMINE,PETSC_DETERMINE,PETSC_DECIDE,dnnz,PETSC_DECIDE,onnz,A,ierr)) 28 PetscCallA(MatSetFromOptions(A,ierr)) 29 PetscCallA(MatSetUp(A,ierr)) 30 deallocate(dnnz) 31 deallocate(onnz) 32 33 !/* This assembly shrinks memory because we do not insert enough number of values */ 34 PetscCallA(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)) 35 PetscCallA(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)) 36 37 !/* MatResetPreallocation restores the memory required by users */ 38 PetscCallA(MatResetPreallocation(A,ierr)) 39 PetscCallA(MatSetOption(A,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_TRUE,ierr)) 40 PetscCallA(MatGetOwnershipRange(A,rstart,rend,ierr)) 41 PetscCallA(MatGetSize(A,M1,N1,ierr)) 42 do i=rstart,rend-1 43 PetscCallA(MatSetValue(A,i,i,two,INSERT_VALUES,ierr)) 44 if (rend<N1) PetscCallA(MatSetValue(A,i,rend,one,INSERT_VALUES,ierr)) 45 end do 46 PetscCallA(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)) 47 PetscCallA(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)) 48 PetscCallA(MatView(A,PETSC_VIEWER_STDOUT_WORLD,ierr)) 49 PetscCallA(MatDestroy(A,ierr)) 50 PetscCallA(PetscFinalize(ierr)) 51 52end program 53 54!/*TEST 55! 56! test: 57! suffix: 1 58! output_file: output/ex4_1.out 59! 60! test: 61! suffix: 2 62! nsize: 2 63! output_file: output/ex4_2.out 64! 65!TEST*/ 66