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