1 /* 2 Formatted test for ISGeneral routines. 3 */ 4 5 static char help[] = "Tests IS general routines.\n\n"; 6 7 #include <petscis.h> 8 #include <petscviewer.h> 9 10 int main(int argc, char **argv) { 11 PetscMPIInt rank, size; 12 PetscInt i, n, *indices; 13 const PetscInt *ii; 14 IS is, newis; 15 PetscBool flg; 16 PetscBool permanent = PETSC_FALSE; 17 PetscBool compute = PETSC_TRUE; 18 19 PetscFunctionBeginUser; 20 PetscCall(PetscInitialize(&argc, &argv, (char *)0, help)); 21 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 22 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 23 24 /* 25 Test IS of size 0 26 */ 27 PetscCall(ISCreateGeneral(PETSC_COMM_SELF, 0, &n, PETSC_COPY_VALUES, &is)); 28 PetscCall(ISGetSize(is, &n)); 29 PetscCheck(n == 0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "ISGetSize"); 30 PetscCall(ISDestroy(&is)); 31 32 /* 33 Create large IS and test ISGetIndices() 34 */ 35 n = 10000 + rank; 36 PetscCall(PetscMalloc1(n, &indices)); 37 for (i = 0; i < n; i++) indices[i] = rank + i; 38 PetscCall(ISCreateGeneral(PETSC_COMM_SELF, n, indices, PETSC_COPY_VALUES, &is)); 39 PetscCall(ISGetIndices(is, &ii)); 40 for (i = 0; i < n; i++) PetscCheck(ii[i] == indices[i], PETSC_COMM_SELF, PETSC_ERR_PLIB, "ISGetIndices"); 41 PetscCall(ISRestoreIndices(is, &ii)); 42 43 /* 44 Check identity and permutation 45 */ 46 /* ISPermutation doesn't check if not set */ 47 PetscCall(ISPermutation(is, &flg)); 48 PetscCheck(!flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "ISPermutation"); 49 PetscCall(ISGetInfo(is, IS_PERMUTATION, IS_LOCAL, compute, &flg)); 50 PetscCheck(rank != 0 || flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "ISGetInfo(IS_PERMUTATION,IS_LOCAL)"); 51 PetscCheck(rank == 0 || !flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "ISGetInfo(IS_PERMUTATION,IS_LOCAL)"); 52 PetscCall(ISIdentity(is, &flg)); 53 PetscCheck(rank != 0 || flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "ISIdentity"); 54 PetscCheck(rank == 0 || !flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "ISIdentity"); 55 PetscCall(ISGetInfo(is, IS_IDENTITY, IS_LOCAL, compute, &flg)); 56 PetscCheck(rank != 0 || flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "ISGetInfo(IS_IDENTITY,IS_LOCAL)"); 57 PetscCheck(rank == 0 || !flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "ISGetInfo(IS_IDENTITY,IS_LOCAL)"); 58 /* we can override the computed values with ISSetInfo() */ 59 PetscCall(ISSetInfo(is, IS_PERMUTATION, IS_LOCAL, permanent, PETSC_TRUE)); 60 PetscCall(ISSetInfo(is, IS_IDENTITY, IS_LOCAL, permanent, PETSC_TRUE)); 61 PetscCall(ISGetInfo(is, IS_PERMUTATION, IS_LOCAL, compute, &flg)); 62 PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "ISGetInfo(IS_PERMUTATION,IS_LOCAL)"); 63 PetscCall(ISGetInfo(is, IS_IDENTITY, IS_LOCAL, compute, &flg)); 64 PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "ISGetInfo(IS_IDENTITY,IS_LOCAL)"); 65 66 PetscCall(ISClearInfoCache(is, PETSC_TRUE)); 67 68 /* 69 Check equality of index sets 70 */ 71 PetscCall(ISEqual(is, is, &flg)); 72 PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "ISEqual"); 73 74 /* 75 Sorting 76 */ 77 PetscCall(ISSort(is)); 78 PetscCall(ISSorted(is, &flg)); 79 PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "ISSort"); 80 PetscCall(ISGetInfo(is, IS_SORTED, IS_LOCAL, compute, &flg)); 81 PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "ISGetInfo(IS_SORTED,IS_LOCAL)"); 82 PetscCall(ISSorted(is, &flg)); 83 PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "ISSort"); 84 PetscCall(ISGetInfo(is, IS_SORTED, IS_LOCAL, compute, &flg)); 85 PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "ISGetInfo(IS_SORTED,IS_LOCAL)"); 86 87 /* 88 Thinks it is a different type? 89 */ 90 PetscCall(PetscObjectTypeCompare((PetscObject)is, ISSTRIDE, &flg)); 91 PetscCheck(!flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "ISStride"); 92 PetscCall(PetscObjectTypeCompare((PetscObject)is, ISBLOCK, &flg)); 93 PetscCheck(!flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "ISBlock"); 94 95 PetscCall(ISDestroy(&is)); 96 97 /* 98 Inverting permutation 99 */ 100 for (i = 0; i < n; i++) indices[i] = n - i - 1; 101 PetscCall(ISCreateGeneral(PETSC_COMM_SELF, n, indices, PETSC_COPY_VALUES, &is)); 102 PetscCall(PetscFree(indices)); 103 PetscCall(ISSetPermutation(is)); 104 PetscCall(ISInvertPermutation(is, PETSC_DECIDE, &newis)); 105 PetscCall(ISGetIndices(newis, &ii)); 106 for (i = 0; i < n; i++) PetscCheck(ii[i] == n - i - 1, PETSC_COMM_SELF, PETSC_ERR_PLIB, "ISInvertPermutation"); 107 PetscCall(ISRestoreIndices(newis, &ii)); 108 PetscCall(ISDestroy(&newis)); 109 PetscCall(ISDestroy(&is)); 110 PetscCall(PetscFinalize()); 111 return 0; 112 } 113 114 /*TEST 115 116 test: 117 nsize: {{1 2 3 4 5}} 118 119 TEST*/ 120