Actual source code: ex1.c
petsc-3.12.1 2019-10-22
1: /*
2: Formatted test for ISGeneral routines.
3: */
5: static char help[] = "Tests IS general routines.\n\n";
7: #include <petscis.h>
8: #include <petscviewer.h>
10: int main(int argc,char **argv)
11: {
12: PetscMPIInt rank,size;
13: PetscInt i,n,*indices;
14: const PetscInt *ii;
15: IS is,newis;
16: PetscBool flg;
19: PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
20: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
21: MPI_Comm_size(PETSC_COMM_WORLD,&size);
23: /*
24: Test IS of size 0
25: */
26: ISCreateGeneral(PETSC_COMM_SELF,0,&n,PETSC_COPY_VALUES,&is);
27: ISGetSize(is,&n);
28: if (n != 0) SETERRQ(PETSC_COMM_SELF,1,"ISGetSize");
29: ISDestroy(&is);
31: /*
32: Create large IS and test ISGetIndices()
33: */
34: n = 10000 + rank;
35: PetscMalloc1(n,&indices);
36: for (i=0; i<n; i++) indices[i] = rank + i;
37: ISCreateGeneral(PETSC_COMM_SELF,n,indices,PETSC_COPY_VALUES,&is);
38: ISGetIndices(is,&ii);
39: for (i=0; i<n; i++) {
40: if (ii[i] != indices[i]) SETERRQ(PETSC_COMM_SELF,1,"ISGetIndices");
41: }
42: ISRestoreIndices(is,&ii);
44: /*
45: Check identity and permutation
46: */
47: ISPermutation(is,&flg);
48: if (flg) SETERRQ(PETSC_COMM_SELF,1,"ISPermutation");
49: ISIdentity(is,&flg);
50: if (!flg) SETERRQ(PETSC_COMM_SELF,1,"ISIdentity");
51: ISSetPermutation(is);
52: ISSetIdentity(is);
53: ISPermutation(is,&flg);
54: if (!flg) SETERRQ(PETSC_COMM_SELF,1,"ISPermutation");
55: ISIdentity(is,&flg);
56: if (!flg) SETERRQ(PETSC_COMM_SELF,1,"ISIdentity");
58: /*
59: Check equality of index sets
60: */
61: ISEqual(is,is,&flg);
62: if (!flg) SETERRQ(PETSC_COMM_SELF,1,"ISEqual");
64: /*
65: Sorting
66: */
67: ISSort(is);
68: ISSorted(is,&flg);
69: if (!flg) SETERRQ(PETSC_COMM_SELF,1,"ISSort");
71: /*
72: Thinks it is a different type?
73: */
74: PetscObjectTypeCompare((PetscObject)is,ISSTRIDE,&flg);
75: if (flg) SETERRQ(PETSC_COMM_SELF,1,"ISStride");
76: PetscObjectTypeCompare((PetscObject)is,ISBLOCK,&flg);
77: if (flg) SETERRQ(PETSC_COMM_SELF,1,"ISBlock");
79: ISDestroy(&is);
81: /*
82: Inverting permutation
83: */
84: for (i=0; i<n; i++) indices[i] = n - i - 1;
85: ISCreateGeneral(PETSC_COMM_SELF,n,indices,PETSC_COPY_VALUES,&is);
86: PetscFree(indices);
87: ISSetPermutation(is);
88: ISInvertPermutation(is,PETSC_DECIDE,&newis);
89: ISGetIndices(newis,&ii);
90: for (i=0; i<n; i++) {
91: if (ii[i] != n - i - 1) SETERRQ(PETSC_COMM_SELF,1,"ISInvertPermutation");
92: }
93: ISRestoreIndices(newis,&ii);
94: ISDestroy(&newis);
95: ISDestroy(&is);
96: PetscFinalize();
97: return ierr;
98: }
100: /*TEST
102: test:
104: TEST*/