1 2 static char help[] = "Tests MatMeshToDual()\n\n"; 3 4 /* 5 Include "petscmat.h" so that we can use matrices. 6 automatically includes: 7 petscsys.h - base PETSc routines petscvec.h - vectors 8 petscmat.h - matrices 9 petscis.h - index sets petscviewer.h - viewers 10 */ 11 #include <petscmat.h> 12 13 int main(int argc,char **args) 14 { 15 Mat mesh,dual; 16 PetscInt Nvertices = 6; /* total number of vertices */ 17 PetscInt ncells = 2; /* number cells on this process */ 18 PetscInt *ii,*jj; 19 PetscMPIInt size,rank; 20 MatPartitioning part; 21 IS is; 22 23 PetscFunctionBeginUser; 24 PetscCall(PetscInitialize(&argc,&args,(char*)0,help)); 25 PetscCallMPI(MPI_Comm_size(MPI_COMM_WORLD,&size)); 26 PetscCheck(size == 2,PETSC_COMM_WORLD,PETSC_ERR_WRONG_MPI_SIZE,"This example is for exactly two processes"); 27 PetscCallMPI(MPI_Comm_rank(MPI_COMM_WORLD,&rank)); 28 29 PetscCall(PetscMalloc1(3,&ii)); 30 PetscCall(PetscMalloc1(6,&jj)); 31 ii[0] = 0; ii[1] = 3; ii[2] = 6; 32 if (rank == 0) { 33 jj[0] = 0; jj[1] = 1; jj[2] = 2; jj[3] = 1; jj[4] = 2; jj[5] = 3; 34 } else { 35 jj[0] = 1; jj[1] = 4; jj[2] = 5; jj[3] = 1; jj[4] = 3; jj[5] = 5; 36 } 37 PetscCall(MatCreateMPIAdj(MPI_COMM_WORLD,ncells,Nvertices,ii,jj,NULL,&mesh)); 38 PetscCall(MatMeshToCellGraph(mesh,2,&dual)); 39 PetscCall(MatView(dual,PETSC_VIEWER_STDOUT_WORLD)); 40 41 PetscCall(MatPartitioningCreate(MPI_COMM_WORLD,&part)); 42 PetscCall(MatPartitioningSetAdjacency(part,dual)); 43 PetscCall(MatPartitioningSetFromOptions(part)); 44 PetscCall(MatPartitioningApply(part,&is)); 45 PetscCall(ISView(is,PETSC_VIEWER_STDOUT_WORLD)); 46 PetscCall(ISDestroy(&is)); 47 PetscCall(MatPartitioningDestroy(&part)); 48 49 PetscCall(MatDestroy(&mesh)); 50 PetscCall(MatDestroy(&dual)); 51 PetscCall(PetscFinalize()); 52 return 0; 53 } 54 55 /*TEST 56 57 build: 58 requires: parmetis 59 60 test: 61 nsize: 2 62 requires: parmetis 63 64 TEST*/ 65