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