1c4762a1bSJed Brown static char help[] = "Example of using graph partitioning to partition a graph\n\n"; 2c4762a1bSJed Brown 3c4762a1bSJed Brown #include <petscmat.h> 4c4762a1bSJed Brown 5c4762a1bSJed Brown int main(int argc, char **args) 6c4762a1bSJed Brown { 7c4762a1bSJed Brown Mat A; 8c4762a1bSJed Brown MatPartitioning part; 9c4762a1bSJed Brown IS is; 10c4762a1bSJed Brown PetscInt r,N = 10, start, end, *vweights; 11c4762a1bSJed Brown PetscBool set_vweights=PETSC_FALSE,use_edge_weights=PETSC_FALSE; 12c4762a1bSJed Brown PetscMPIInt rank; 13c4762a1bSJed Brown MPI_Comm comm; 14c4762a1bSJed Brown 15*327415f7SBarry Smith PetscFunctionBeginUser; 169566063dSJacob Faibussowitsch PetscCall(PetscInitialize(&argc, &args, (char*) 0, help)); 17c4762a1bSJed Brown comm = PETSC_COMM_WORLD; 189566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetInt(NULL,NULL, "-N", &N, NULL)); 199566063dSJacob Faibussowitsch PetscCallMPI(MPI_Comm_rank(comm,&rank)); 209566063dSJacob Faibussowitsch PetscCall(MatCreate(comm, &A)); 219566063dSJacob Faibussowitsch PetscCall(MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, N, N)); 229566063dSJacob Faibussowitsch PetscCall(MatSetFromOptions(A)); 239566063dSJacob Faibussowitsch PetscCall(MatSeqAIJSetPreallocation(A, 3, NULL)); 249566063dSJacob Faibussowitsch PetscCall(MatMPIAIJSetPreallocation(A, 3, NULL, 2, NULL)); 259566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetBool(NULL,NULL,"-test_vertex_weights",&set_vweights,NULL)); 269566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetBool(NULL,NULL,"-test_use_edge_weights",&use_edge_weights,NULL)); 27c4762a1bSJed Brown /* Create a linear mesh */ 289566063dSJacob Faibussowitsch PetscCall(MatGetOwnershipRange(A, &start, &end)); 29c4762a1bSJed Brown if (set_vweights) { 309566063dSJacob Faibussowitsch PetscCall(PetscMalloc1(end-start,&vweights)); 31c4762a1bSJed Brown for (r = start; r < end; ++r) 32c4762a1bSJed Brown vweights[r-start] = rank+1; 33c4762a1bSJed Brown } 34c4762a1bSJed Brown for (r = start; r < end; ++r) { 35c4762a1bSJed Brown if (r == 0) { 36c4762a1bSJed Brown PetscInt cols[2]; 37c4762a1bSJed Brown PetscScalar vals[2]; 38c4762a1bSJed Brown 39c4762a1bSJed Brown cols[0] = r; cols[1] = r+1; 40c4762a1bSJed Brown vals[0] = 1.0; vals[1] = use_edge_weights? 2.0: 1.0; 41c4762a1bSJed Brown 429566063dSJacob Faibussowitsch PetscCall(MatSetValues(A, 1, &r, 2, cols, vals, INSERT_VALUES)); 43c4762a1bSJed Brown } else if (r == N-1) { 44c4762a1bSJed Brown PetscInt cols[2]; 45c4762a1bSJed Brown PetscScalar vals[2]; 46c4762a1bSJed Brown 47c4762a1bSJed Brown cols[0] = r-1; cols[1] = r; 48c4762a1bSJed Brown vals[0] = use_edge_weights? 3.0:1.0; vals[1] = 1.0; 49c4762a1bSJed Brown 509566063dSJacob Faibussowitsch PetscCall(MatSetValues(A, 1, &r, 2, cols, vals, INSERT_VALUES)); 51c4762a1bSJed Brown } else { 52c4762a1bSJed Brown PetscInt cols[3]; 53c4762a1bSJed Brown PetscScalar vals[3]; 54c4762a1bSJed Brown 55c4762a1bSJed Brown cols[0] = r-1; cols[1] = r; cols[2] = r+1; 56c4762a1bSJed Brown /* ADJ matrix needs to be symmetric */ 57c4762a1bSJed Brown vals[0] = use_edge_weights? (cols[0]==0? 2.0:5.0):1.0; 58c4762a1bSJed Brown vals[1] = 1.0; 59c4762a1bSJed Brown vals[2] = use_edge_weights? (cols[2]==N-1? 3.0:5.0):1.0; 60c4762a1bSJed Brown 619566063dSJacob Faibussowitsch PetscCall(MatSetValues(A, 1, &r, 3, cols, vals, INSERT_VALUES)); 62c4762a1bSJed Brown } 63c4762a1bSJed Brown } 649566063dSJacob Faibussowitsch PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 659566063dSJacob Faibussowitsch PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 66c4762a1bSJed Brown 679566063dSJacob Faibussowitsch PetscCall(MatPartitioningCreate(comm, &part)); 689566063dSJacob Faibussowitsch PetscCall(MatPartitioningSetAdjacency(part, A)); 691baa6e33SBarry Smith if (set_vweights) PetscCall(MatPartitioningSetVertexWeights(part,vweights)); 70c4762a1bSJed Brown if (use_edge_weights) { 719566063dSJacob Faibussowitsch PetscCall(MatPartitioningSetUseEdgeWeights(part,use_edge_weights)); 72c4762a1bSJed Brown 739566063dSJacob Faibussowitsch PetscCall(MatPartitioningGetUseEdgeWeights(part,&use_edge_weights)); 7428b400f6SJacob Faibussowitsch PetscCheck(use_edge_weights,comm,PETSC_ERR_ARG_INCOMP, "use_edge_weights flag does not setup correctly "); 75c4762a1bSJed Brown } 769566063dSJacob Faibussowitsch PetscCall(MatPartitioningSetFromOptions(part)); 779566063dSJacob Faibussowitsch PetscCall(MatPartitioningApply(part, &is)); 789566063dSJacob Faibussowitsch PetscCall(ISView(is, PETSC_VIEWER_STDOUT_WORLD)); 799566063dSJacob Faibussowitsch PetscCall(ISDestroy(&is)); 809566063dSJacob Faibussowitsch PetscCall(MatPartitioningDestroy(&part)); 81c4762a1bSJed Brown 829566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A)); 839566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 84b122ec5aSJacob Faibussowitsch return 0; 85c4762a1bSJed Brown } 86c4762a1bSJed Brown 87c4762a1bSJed Brown /*TEST 88c4762a1bSJed Brown 89c4762a1bSJed Brown test: 90c4762a1bSJed Brown nsize: 3 91c4762a1bSJed Brown requires: parmetis 92c4762a1bSJed Brown args: -mat_partitioning_type parmetis 93c4762a1bSJed Brown 94c4762a1bSJed Brown test: 95c4762a1bSJed Brown suffix: 2 96c4762a1bSJed Brown nsize: 3 97c4762a1bSJed Brown requires: ptscotch 98c4762a1bSJed Brown args: -mat_partitioning_type ptscotch 99c4762a1bSJed Brown 100c4762a1bSJed Brown test: 101c4762a1bSJed Brown suffix: 3 102c4762a1bSJed Brown nsize: 4 103c4762a1bSJed Brown requires: party 104c4762a1bSJed Brown args: -mat_partitioning_type party 105c4762a1bSJed Brown 106c4762a1bSJed Brown test: 107c4762a1bSJed Brown suffix: 4 108c4762a1bSJed Brown nsize: 3 109c4762a1bSJed Brown requires: chaco 110c4762a1bSJed Brown args: -mat_partitioning_type chaco 111c4762a1bSJed Brown 112c4762a1bSJed Brown test: 113c4762a1bSJed Brown suffix: 5 114c4762a1bSJed Brown nsize: 3 115c4762a1bSJed Brown requires: parmetis 116c4762a1bSJed Brown args: -mat_partitioning_type hierarch -mat_partitioning_hierarchical_nfineparts 3 -mat_partitioning_nparts 10 -N 100 117c4762a1bSJed Brown 118c4762a1bSJed Brown test: 119c4762a1bSJed Brown suffix: 6 120c4762a1bSJed Brown nsize: 3 121c4762a1bSJed Brown requires: parmetis 122c4762a1bSJed Brown args: -mat_partitioning_type hierarch -mat_partitioning_hierarchical_nfineparts 3 -mat_partitioning_nparts 10 -N 100 -test_vertex_weights 1 -mat_partitioning_use_edge_weights 1 123c4762a1bSJed Brown 124c4762a1bSJed Brown test: 125c4762a1bSJed Brown suffix: 7 126c4762a1bSJed Brown nsize: 2 127c4762a1bSJed Brown requires: parmetis 128c4762a1bSJed Brown args: -mat_partitioning_type hierarch -mat_partitioning_hierarchical_nfineparts 2 -mat_partitioning_nparts 10 -mat_partitioning_hierarchical_fineparttype hierarch -malloc_dump -N 100 -mat_partitioning_improve 1 129c4762a1bSJed Brown 130c4762a1bSJed Brown test: 131c4762a1bSJed Brown suffix: 8 132c4762a1bSJed Brown nsize: 2 133c4762a1bSJed Brown requires: parmetis 134c4762a1bSJed Brown args: -mat_partitioning_type parmetis -mat_partitioning_nparts 3 -test_use_edge_weights 1 135c4762a1bSJed Brown 136c4762a1bSJed Brown test: 137c4762a1bSJed Brown suffix: 9 138c4762a1bSJed Brown nsize: 2 139c4762a1bSJed Brown requires: ptscotch 140c4762a1bSJed Brown args: -mat_partitioning_type ptscotch -mat_partitioning_nparts 3 -test_use_edge_weights 1 -mat_partitioning_ptscotch_proc_weight 0 141c4762a1bSJed Brown 142c4762a1bSJed Brown TEST*/ 143