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