1*c4762a1bSJed Brown /* 2*c4762a1bSJed Brown * ex193.c 3*c4762a1bSJed Brown * 4*c4762a1bSJed Brown * Created on: Jul 29, 2015 5*c4762a1bSJed Brown * Author: Fande Kong fdkong.jd@gmail.com 6*c4762a1bSJed Brown */ 7*c4762a1bSJed Brown /* 8*c4762a1bSJed Brown * An example demonstrates how to use hierarchical partitioning approach 9*c4762a1bSJed Brown */ 10*c4762a1bSJed Brown 11*c4762a1bSJed Brown #include <petscmat.h> 12*c4762a1bSJed Brown 13*c4762a1bSJed Brown static char help[] = "Illustrates use of hierarchical partitioning.\n"; 14*c4762a1bSJed Brown 15*c4762a1bSJed Brown int main(int argc,char **args) 16*c4762a1bSJed Brown { 17*c4762a1bSJed Brown Mat A; /* matrix */ 18*c4762a1bSJed Brown PetscInt m,n; /* mesh dimensions in x- and y- directions */ 19*c4762a1bSJed Brown PetscInt i,j,Ii,J,Istart,Iend; 20*c4762a1bSJed Brown PetscErrorCode ierr; 21*c4762a1bSJed Brown PetscMPIInt size; 22*c4762a1bSJed Brown PetscScalar v; 23*c4762a1bSJed Brown MatPartitioning part; 24*c4762a1bSJed Brown IS coarseparts,fineparts; 25*c4762a1bSJed Brown IS is,isn,isrows; 26*c4762a1bSJed Brown MPI_Comm comm; 27*c4762a1bSJed Brown 28*c4762a1bSJed Brown ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr; 29*c4762a1bSJed Brown comm = PETSC_COMM_WORLD; 30*c4762a1bSJed Brown ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr); 31*c4762a1bSJed Brown ierr = PetscOptionsBegin(comm,NULL,"ex193","hierarchical partitioning");CHKERRQ(ierr); 32*c4762a1bSJed Brown m = 15; 33*c4762a1bSJed Brown ierr = PetscOptionsInt("-M","Number of mesh points in the x-direction","partitioning",m,&m,NULL);CHKERRQ(ierr); 34*c4762a1bSJed Brown n = 15; 35*c4762a1bSJed Brown ierr = PetscOptionsInt("-N","Number of mesh points in the y-direction","partitioning",n,&n,NULL);CHKERRQ(ierr); 36*c4762a1bSJed Brown ierr = PetscOptionsEnd();CHKERRQ(ierr); 37*c4762a1bSJed Brown 38*c4762a1bSJed Brown /* 39*c4762a1bSJed Brown Assemble the matrix for the five point stencil (finite difference), YET AGAIN 40*c4762a1bSJed Brown */ 41*c4762a1bSJed Brown ierr = MatCreate(comm,&A);CHKERRQ(ierr); 42*c4762a1bSJed Brown ierr = MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,m*n,m*n);CHKERRQ(ierr); 43*c4762a1bSJed Brown ierr = MatSetFromOptions(A);CHKERRQ(ierr); 44*c4762a1bSJed Brown ierr = MatSetUp(A);CHKERRQ(ierr); 45*c4762a1bSJed Brown ierr = MatGetOwnershipRange(A,&Istart,&Iend);CHKERRQ(ierr); 46*c4762a1bSJed Brown for (Ii=Istart; Ii<Iend; Ii++) { 47*c4762a1bSJed Brown v = -1.0; i = Ii/n; j = Ii - i*n; 48*c4762a1bSJed Brown if (i>0) {J = Ii - n; ierr = MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);CHKERRQ(ierr);} 49*c4762a1bSJed Brown if (i<m-1) {J = Ii + n; ierr = MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);CHKERRQ(ierr);} 50*c4762a1bSJed Brown if (j>0) {J = Ii - 1; ierr = MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);CHKERRQ(ierr);} 51*c4762a1bSJed Brown if (j<n-1) {J = Ii + 1; ierr = MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);CHKERRQ(ierr);} 52*c4762a1bSJed Brown v = 4.0; ierr = MatSetValues(A,1,&Ii,1,&Ii,&v,INSERT_VALUES);CHKERRQ(ierr); 53*c4762a1bSJed Brown } 54*c4762a1bSJed Brown ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 55*c4762a1bSJed Brown ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 56*c4762a1bSJed Brown ierr = MatView(A,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 57*c4762a1bSJed Brown /* 58*c4762a1bSJed Brown Partition the graph of the matrix 59*c4762a1bSJed Brown */ 60*c4762a1bSJed Brown ierr = MatPartitioningCreate(comm,&part);CHKERRQ(ierr); 61*c4762a1bSJed Brown ierr = MatPartitioningSetAdjacency(part,A);CHKERRQ(ierr); 62*c4762a1bSJed Brown ierr = MatPartitioningSetType(part,MATPARTITIONINGHIERARCH);CHKERRQ(ierr); 63*c4762a1bSJed Brown ierr = MatPartitioningHierarchicalSetNcoarseparts(part,2);CHKERRQ(ierr); 64*c4762a1bSJed Brown ierr = MatPartitioningHierarchicalSetNfineparts(part,4);CHKERRQ(ierr); 65*c4762a1bSJed Brown ierr = MatPartitioningSetFromOptions(part);CHKERRQ(ierr); 66*c4762a1bSJed Brown /* get new processor owner number of each vertex */ 67*c4762a1bSJed Brown ierr = MatPartitioningApply(part,&is);CHKERRQ(ierr); 68*c4762a1bSJed Brown /* coarse parts */ 69*c4762a1bSJed Brown ierr = MatPartitioningHierarchicalGetCoarseparts(part,&coarseparts);CHKERRQ(ierr); 70*c4762a1bSJed Brown ierr = ISView(coarseparts,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 71*c4762a1bSJed Brown /* fine parts */ 72*c4762a1bSJed Brown ierr = MatPartitioningHierarchicalGetFineparts(part,&fineparts);CHKERRQ(ierr); 73*c4762a1bSJed Brown ierr = ISView(fineparts,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 74*c4762a1bSJed Brown /* partitioning */ 75*c4762a1bSJed Brown ierr = ISView(is,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 76*c4762a1bSJed Brown /* get new global number of each old global number */ 77*c4762a1bSJed Brown ierr = ISPartitioningToNumbering(is,&isn);CHKERRQ(ierr); 78*c4762a1bSJed Brown ierr = ISView(isn,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 79*c4762a1bSJed Brown ierr = ISBuildTwoSided(is,NULL,&isrows);CHKERRQ(ierr); 80*c4762a1bSJed Brown ierr = ISView(isrows,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 81*c4762a1bSJed Brown ierr = ISDestroy(&is);CHKERRQ(ierr); 82*c4762a1bSJed Brown ierr = ISDestroy(&coarseparts);CHKERRQ(ierr); 83*c4762a1bSJed Brown ierr = ISDestroy(&fineparts);CHKERRQ(ierr); 84*c4762a1bSJed Brown ierr = ISDestroy(&isrows);CHKERRQ(ierr); 85*c4762a1bSJed Brown ierr = ISDestroy(&isn);CHKERRQ(ierr); 86*c4762a1bSJed Brown ierr = MatPartitioningDestroy(&part);CHKERRQ(ierr); 87*c4762a1bSJed Brown ierr = MatDestroy(&A);CHKERRQ(ierr); 88*c4762a1bSJed Brown ierr = PetscFinalize(); 89*c4762a1bSJed Brown return ierr; 90*c4762a1bSJed Brown } 91*c4762a1bSJed Brown 92*c4762a1bSJed Brown 93*c4762a1bSJed Brown 94*c4762a1bSJed Brown 95*c4762a1bSJed Brown /*TEST 96*c4762a1bSJed Brown 97*c4762a1bSJed Brown test: 98*c4762a1bSJed Brown nsize: 4 99*c4762a1bSJed Brown args: -mat_partitioning_hierarchical_Nfineparts 2 100*c4762a1bSJed Brown requires: parmetis 101*c4762a1bSJed Brown TODO: cannot run because parmetis does reproduce across all machines, probably due to nonportable random number generator 102*c4762a1bSJed Brown 103*c4762a1bSJed Brown TEST*/ 104