1 static char help[] = "Example of using graph partitioning with a matrix in which some procs have empty ownership\n\n"; 2 3 /*T 4 Concepts: Mat^mat partitioning 5 Concepts: Mat^image segmentation 6 Processors: n 7 T*/ 8 9 #include <petscmat.h> 10 11 int main(int argc, char **args) 12 { 13 Mat A; 14 MatPartitioning part; 15 IS is; 16 PetscInt i,m,N,rstart,rend,nemptyranks,*emptyranks,nbigranks,*bigranks; 17 PetscMPIInt rank,size; 18 PetscErrorCode ierr; 19 20 PetscCall(PetscInitialize(&argc,&args,(char*)0,help)); 21 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 22 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 23 24 nemptyranks = 10; 25 nbigranks = 10; 26 PetscCall(PetscMalloc2(nemptyranks,&emptyranks,nbigranks,&bigranks)); 27 28 ierr = PetscOptionsBegin(PETSC_COMM_WORLD,NULL,"Partitioning example options",NULL);PetscCall(ierr); 29 PetscCall(PetscOptionsIntArray("-emptyranks","Ranks to be skipped by partition","",emptyranks,&nemptyranks,NULL)); 30 PetscCall(PetscOptionsIntArray("-bigranks","Ranks to be overloaded","",bigranks,&nbigranks,NULL)); 31 ierr = PetscOptionsEnd();PetscCall(ierr); 32 33 m = 1; 34 for (i=0; i<nemptyranks; i++) if (rank == emptyranks[i]) m = 0; 35 for (i=0; i<nbigranks; i++) if (rank == bigranks[i]) m = 5; 36 37 PetscCall(MatCreate(PETSC_COMM_WORLD,&A)); 38 PetscCall(MatSetSizes(A,m,m,PETSC_DECIDE,PETSC_DECIDE)); 39 PetscCall(MatSetFromOptions(A)); 40 PetscCall(MatSeqAIJSetPreallocation(A,3,NULL)); 41 PetscCall(MatMPIAIJSetPreallocation(A,3,NULL,2,NULL)); 42 PetscCall(MatSeqBAIJSetPreallocation(A,1,3,NULL)); 43 PetscCall(MatMPIBAIJSetPreallocation(A,1,3,NULL,2,NULL)); 44 PetscCall(MatSeqSBAIJSetPreallocation(A,1,2,NULL)); 45 PetscCall(MatMPISBAIJSetPreallocation(A,1,2,NULL,1,NULL)); 46 47 PetscCall(MatGetSize(A,NULL,&N)); 48 PetscCall(MatGetOwnershipRange(A,&rstart,&rend)); 49 for (i=rstart; i<rend; i++) { 50 const PetscInt cols[] = {(i+N-1)%N,i,(i+1)%N}; 51 const PetscScalar vals[] = {1,1,1}; 52 PetscCall(MatSetValues(A,1,&i,3,cols,vals,INSERT_VALUES)); 53 } 54 PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY)); 55 PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY)); 56 PetscCall(MatView(A,PETSC_VIEWER_STDOUT_WORLD)); 57 58 PetscCall(MatPartitioningCreate(PETSC_COMM_WORLD,&part)); 59 PetscCall(MatPartitioningSetAdjacency(part,A)); 60 PetscCall(MatPartitioningSetFromOptions(part)); 61 PetscCall(MatPartitioningApply(part,&is)); 62 PetscCall(ISView(is,PETSC_VIEWER_STDOUT_WORLD)); 63 PetscCall(ISDestroy(&is)); 64 PetscCall(MatPartitioningDestroy(&part)); 65 PetscCall(MatDestroy(&A)); 66 PetscCall(PetscFree2(emptyranks,bigranks)); 67 PetscCall(PetscFinalize()); 68 return 0; 69 } 70 71 /*TEST 72 73 test: 74 nsize: 8 75 args: -emptyranks 0,2,4 -bigranks 1,3,7 -mat_partitioning_type average 76 # cannot test with external package partitioners since they produce different results on different systems 77 78 TEST*/ 79