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