1 static char help[] = "Tests ISRenumber.\n\n"; 2 3 #include <petscis.h> 4 5 PetscErrorCode TestRenumber(IS is, IS mult) 6 { 7 IS nis; 8 PetscInt N; 9 10 PetscFunctionBegin; 11 CHKERRQ(PetscPrintf(PetscObjectComm((PetscObject)is),"\n-----------------\n")); 12 CHKERRQ(PetscPrintf(PetscObjectComm((PetscObject)is),"\nInitial\n")); 13 CHKERRQ(ISView(is,NULL)); 14 if (mult) { 15 CHKERRQ(PetscPrintf(PetscObjectComm((PetscObject)is),"\nMult\n")); 16 CHKERRQ(ISView(mult,NULL)); 17 } 18 CHKERRQ(ISRenumber(is,mult,&N,NULL)); 19 CHKERRQ(PetscPrintf(PetscObjectComm((PetscObject)is),"\nRenumbered, unique entries %" PetscInt_FMT "\n",N)); 20 CHKERRQ(ISRenumber(is,mult,NULL,&nis)); 21 CHKERRQ(ISView(nis,NULL)); 22 CHKERRQ(ISDestroy(&nis)); 23 PetscFunctionReturn(0); 24 } 25 26 int main(int argc, char **argv) 27 { 28 IS is; 29 PetscMPIInt size, rank; 30 31 CHKERRQ(PetscInitialize(&argc, &argv, NULL, help)); 32 CHKERRMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 33 CHKERRMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 34 35 for (PetscInt c = 0; c < 3; c++) { 36 IS mult = NULL; 37 38 CHKERRQ(ISCreateStride(PETSC_COMM_WORLD,0,0,0,&is)); 39 if (c) { 40 PetscInt n; 41 CHKERRQ(ISGetLocalSize(is,&n)); 42 CHKERRQ(ISCreateStride(PETSC_COMM_WORLD,n,c-2,0,&mult)); 43 } 44 CHKERRQ(TestRenumber(is,mult)); 45 CHKERRQ(ISDestroy(&is)); 46 CHKERRQ(ISDestroy(&mult)); 47 48 CHKERRQ(ISCreateStride(PETSC_COMM_WORLD,2,-rank-1,-4,&is)); 49 if (c) { 50 PetscInt n; 51 CHKERRQ(ISGetLocalSize(is,&n)); 52 CHKERRQ(ISCreateStride(PETSC_COMM_WORLD,n,c-2,0,&mult)); 53 } 54 CHKERRQ(TestRenumber(is,mult)); 55 CHKERRQ(ISDestroy(&is)); 56 CHKERRQ(ISDestroy(&mult)); 57 58 CHKERRQ(ISCreateStride(PETSC_COMM_WORLD,10,4+rank,2,&is)); 59 if (c) { 60 PetscInt n; 61 CHKERRQ(ISGetLocalSize(is,&n)); 62 CHKERRQ(ISCreateStride(PETSC_COMM_WORLD,n,c-2,1,&mult)); 63 } 64 CHKERRQ(TestRenumber(is,mult)); 65 CHKERRQ(ISDestroy(&is)); 66 CHKERRQ(ISDestroy(&mult)); 67 68 CHKERRQ(ISCreateStride(PETSC_COMM_WORLD,10,-rank-1,2,&is)); 69 if (c) { 70 PetscInt n; 71 CHKERRQ(ISGetLocalSize(is,&n)); 72 CHKERRQ(ISCreateStride(PETSC_COMM_WORLD,n,c-2,1,&mult)); 73 } 74 CHKERRQ(TestRenumber(is,mult)); 75 CHKERRQ(ISDestroy(&is)); 76 CHKERRQ(ISDestroy(&mult)); 77 } 78 /* Finalize */ 79 CHKERRQ(PetscFinalize()); 80 return 0; 81 } 82 83 /*TEST 84 85 test: 86 suffix: 1 87 nsize: {{1 2}separate output} 88 89 TEST*/ 90