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