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