1*c4762a1bSJed Brown 2*c4762a1bSJed Brownprogram main 3*c4762a1bSJed Brown#include <petsc/finclude/petscvec.h> 4*c4762a1bSJed Brown#include <petsc/finclude/petscmat.h> 5*c4762a1bSJed Brown 6*c4762a1bSJed Brownuse petscvec 7*c4762a1bSJed Brownuse petscmat 8*c4762a1bSJed Brown 9*c4762a1bSJed Brownimplicit none 10*c4762a1bSJed Brown 11*c4762a1bSJed Brown Mat A 12*c4762a1bSJed Brown MatPartitioning part 13*c4762a1bSJed Brown IS is 14*c4762a1bSJed Brown PetscInt :: i,m,N 15*c4762a1bSJed Brown PetscInt :: rstart,rend 16*c4762a1bSJed Brown PetscInt,pointer,dimension(:) :: emptyranks,bigranks,cols 17*c4762a1bSJed Brown PetscScalar,pointer,dimension(:) :: vals 18*c4762a1bSJed Brown PetscInt :: & 19*c4762a1bSJed Brown nbigranks = 10, & 20*c4762a1bSJed Brown nemptyranks = 10 21*c4762a1bSJed Brown PetscMPIInt :: rank,sizef 22*c4762a1bSJed Brown PetscErrorCode ierr 23*c4762a1bSJed Brown PetscBool set 24*c4762a1bSJed Brown PetscInt,parameter :: zero = 0, one = 1, two = 2, three = 3 25*c4762a1bSJed Brown 26*c4762a1bSJed Brown call PetscInitialize(PETSC_NULL_CHARACTER,ierr) 27*c4762a1bSJed Brown if (ierr /= 0) then 28*c4762a1bSJed Brown print*,'PetscInitialize failed' 29*c4762a1bSJed Brown stop 30*c4762a1bSJed Brown endif 31*c4762a1bSJed Brown 32*c4762a1bSJed Brown call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr);CHKERRA(ierr) 33*c4762a1bSJed Brown call MPI_Comm_size(PETSC_COMM_WORLD,sizef,ierr);CHKERRA(ierr) 34*c4762a1bSJed Brown 35*c4762a1bSJed Brown allocate(emptyranks(nemptyranks)) 36*c4762a1bSJed Brown allocate(bigranks(nbigranks)) 37*c4762a1bSJed Brown 38*c4762a1bSJed Brown call PetscOptionsGetIntArray(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,"-emptyranks",emptyranks,nemptyranks,set,ierr);CHKERRA(ierr) 39*c4762a1bSJed Brown call PetscOptionsGetIntArray(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,"-bigranks",bigranks,nbigranks,set,ierr);CHKERRA(ierr) 40*c4762a1bSJed Brown 41*c4762a1bSJed Brown m = 1 42*c4762a1bSJed Brown do i=1,nemptyranks 43*c4762a1bSJed Brown if (rank == emptyranks(i)) m = 0 44*c4762a1bSJed Brown end do 45*c4762a1bSJed Brown do i=1,nbigranks 46*c4762a1bSJed Brown if (rank == bigranks(i)) m = 5 47*c4762a1bSJed Brown end do 48*c4762a1bSJed Brown 49*c4762a1bSJed Brown deallocate(emptyranks) 50*c4762a1bSJed Brown deallocate(bigranks) 51*c4762a1bSJed Brown 52*c4762a1bSJed Brown call MatCreate(PETSC_COMM_WORLD,A,ierr);CHKERRA(ierr) 53*c4762a1bSJed Brown call MatSetsizes(A,m,m,PETSC_DECIDE,PETSC_DECIDE,ierr);CHKERRA(ierr) 54*c4762a1bSJed Brown call MatSetFromOptions(A,ierr);CHKERRA(ierr) 55*c4762a1bSJed Brown call MatSeqAIJSetPreallocation(A,three,PETSC_NULL_INTEGER,ierr);CHKERRA(ierr) 56*c4762a1bSJed Brown call MatMPIAIJSetPreallocation(A,three,PETSC_NULL_INTEGER,two,PETSC_NULL_INTEGER,ierr);CHKERRA(ierr) 57*c4762a1bSJed Brown call MatSeqBAIJSetPreallocation(A,one,three,PETSC_NULL_INTEGER,ierr);CHKERRA(ierr) 58*c4762a1bSJed Brown call MatMPIBAIJSetPreallocation(A,one,three,PETSC_NULL_INTEGER,2,PETSC_NULL_INTEGER,ierr);CHKERRA(ierr) 59*c4762a1bSJed Brown call MatSeqSBAIJSetPreallocation(A,one,two,PETSC_NULL_INTEGER,ierr);CHKERRA(ierr) 60*c4762a1bSJed Brown call MatMPISBAIJSetPreallocation(A,one,two,PETSC_NULL_INTEGER,1,PETSC_NULL_INTEGER,ierr);CHKERRA(ierr) 61*c4762a1bSJed Brown 62*c4762a1bSJed Brown call MatGetSize(A,PETSC_NULL_INTEGER,N,ierr);CHKERRA(ierr) 63*c4762a1bSJed Brown call MatGetOwnershipRange(A,rstart,rend,ierr);CHKERRA(ierr) 64*c4762a1bSJed Brown 65*c4762a1bSJed Brown allocate(cols(0:3)) 66*c4762a1bSJed Brown allocate(vals(0:3)) 67*c4762a1bSJed Brown do i=rstart,rend-1 68*c4762a1bSJed Brown 69*c4762a1bSJed Brown cols = (/mod((i+N-1),N),i,mod((i+1),N)/) 70*c4762a1bSJed Brown vals = [1.0,1.0,1.0] 71*c4762a1bSJed Brown call MatSetValues(A,one,i,three,cols,vals,INSERT_VALUES,ierr);CHKERRA(ierr) 72*c4762a1bSJed Brown end do 73*c4762a1bSJed Brown deallocate(cols) 74*c4762a1bSJed Brown deallocate(vals) 75*c4762a1bSJed Brown call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr);CHKERRA(ierr) 76*c4762a1bSJed Brown call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr);CHKERRA(ierr) 77*c4762a1bSJed Brown call MatView(A,PETSC_VIEWER_STDOUT_WORLD,ierr);CHKERRA(ierr) 78*c4762a1bSJed Brown 79*c4762a1bSJed Brown call MatPartitioningCreate(PETSC_COMM_WORLD,part,ierr);CHKERRA(ierr) 80*c4762a1bSJed Brown call MatPartitioningSetAdjacency(part,A,ierr);CHKERRA(ierr) 81*c4762a1bSJed Brown call MatPartitioningSetFromOptions(part,ierr);CHKERRA(ierr) 82*c4762a1bSJed Brown call MatPartitioningApply(part,is,ierr);CHKERRA(ierr) 83*c4762a1bSJed Brown call ISView(is,PETSC_VIEWER_STDOUT_WORLD,ierr);CHKERRA(ierr) 84*c4762a1bSJed Brown call ISDestroy(is,ierr);CHKERRA(ierr) 85*c4762a1bSJed Brown call MatPartitioningDestroy(part,ierr);CHKERRA(ierr) 86*c4762a1bSJed Brown call MatDestroy(A,ierr);CHKERRA(ierr) 87*c4762a1bSJed Brown call PetscFinalize(ierr);CHKERRA(ierr) 88*c4762a1bSJed Brown 89*c4762a1bSJed Brownend program 90*c4762a1bSJed Brown 91*c4762a1bSJed Brown 92*c4762a1bSJed Brown!/*TEST 93*c4762a1bSJed Brown! 94*c4762a1bSJed Brown! test: 95*c4762a1bSJed Brown! nsize: 8 96*c4762a1bSJed Brown! args: -emptyranks 0,2,4 -bigranks 1,3,7 -mat_partitioning_type average 97*c4762a1bSJed Brown! output_file: output/ex17_1.out 98*c4762a1bSJed Brown! # cannot test with external package partitioners since they produce different results on different systems 99*c4762a1bSJed Brown! 100*c4762a1bSJed Brown!TEST*/ 101