1program main 2#include <petsc/finclude/petscvec.h> 3#include <petsc/finclude/petscmat.h> 4 5 use petscvec 6 use petscmat 7 8 implicit 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