1#include <petsc/finclude/petscmat.h> 2program main 3 use petscmat 4 5 implicit none 6 7 Mat :: A 8 MatPartitioning :: part 9 IS :: is 10 PetscInt :: r, myStart, myEnd 11 PetscInt :: N = 10 12 PetscErrorCode :: ierr 13 PetscScalar, pointer, dimension(:) :: vals 14 PetscInt, pointer, dimension(:) :: cols 15 PetscBool :: flg 16 PetscInt, parameter :: one = 1, two = 2, three = 3 17 18 PetscCallA(PetscInitialize(ierr)) 19 20 PetscCallA(PetscOptionsGetInt(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-N', N, flg, ierr)) 21 PetscCallA(MatCreate(PETSC_COMM_WORLD, A, ierr)) 22 PetscCallA(MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, N, N, ierr)) 23 PetscCallA(MatSetFromOptions(A, ierr)) 24 PetscCallA(MatSeqAIJSetPreallocation(A, three, PETSC_NULL_INTEGER_ARRAY, ierr)) 25 PetscCallA(MatMPIAIJSetPreallocation(A, three, PETSC_NULL_INTEGER_ARRAY, two, PETSC_NULL_INTEGER_ARRAY, ierr)) 26 27 !/* Create a linear mesh */ 28 PetscCallA(MatGetOwnershipRange(A, myStart, myEnd, ierr)) 29 30 do r = myStart, myEnd - 1 31 if (r == 0) then 32 allocate (vals(2)) 33 vals = 1.0 34 allocate (cols(2), source=[r, r + 1]) 35 PetscCallA(MatSetValues(A, one, [r], two, cols, vals, INSERT_VALUES, ierr)) 36 deallocate (cols) 37 deallocate (vals) 38 else if (r == N - 1) then 39 allocate (vals(2)) 40 vals = 1.0 41 allocate (cols(2), source=[r - 1, r]) 42 PetscCallA(MatSetValues(A, one, [r], two, cols, vals, INSERT_VALUES, ierr)) 43 deallocate (cols) 44 deallocate (vals) 45 else 46 allocate (vals(3)) 47 vals = 1.0 48 allocate (cols(3), source=[r - 1, r, r + 1]) 49 PetscCallA(MatSetValues(A, one, [r], three, cols, vals, INSERT_VALUES, ierr)) 50 deallocate (cols) 51 deallocate (vals) 52 end if 53 end do 54 PetscCallA(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY, ierr)) 55 PetscCallA(MatAssemblyend(A, MAT_FINAL_ASSEMBLY, ierr)) 56 PetscCallA(MatPartitioningCreate(PETSC_COMM_WORLD, part, ierr)) 57 PetscCallA(MatPartitioningSetAdjacency(part, A, ierr)) 58 PetscCallA(MatPartitioningSetFromOptions(part, ierr)) 59 PetscCallA(MatPartitioningApply(part, is, ierr)) 60 PetscCallA(ISView(is, PETSC_VIEWER_STDOUT_WORLD, ierr)) 61 PetscCallA(ISDestroy(is, ierr)) 62 PetscCallA(MatPartitioningDestroy(part, ierr)) 63 PetscCallA(MatDestroy(A, ierr)) 64 PetscCallA(PetscFinalize(ierr)) 65 66end program 67 68!/*TEST 69! 70! test: 71! nsize: 3 72! requires: parmetis 73! args: -mat_partitioning_type parmetis 74! output_file: output/ex15_1.out 75! 76! test: 77! suffix: 2 78! nsize: 3 79! requires: ptscotch 80! args: -mat_partitioning_type ptscotch 81! output_file: output/ex15_2.out 82! 83! test: 84! suffix: 3 85! nsize: 4 86! requires: party 87! args: -mat_partitioning_type party 88! output_file: output/ex15_3.out 89! 90! test: 91! suffix: 4 92! nsize: 3 93! requires: chaco 94! args: -mat_partitioning_type chaco 95! output_file: output/ex15_4.out 96! 97! test: 98! suffix: 5 99! nsize: 3 100! requires: parmetis 101! args: -mat_partitioning_type hierarch -mat_partitioning_hierarchical_nfineparts 3 -mat_partitioning_nparts 10 -N 100 102! output_file: output/ex15_5.out 103! 104!TEST*/ 105