xref: /petsc/src/sys/tests/ex40.c (revision 1e1ea65d8de51fde77ce8a787efbef25e407badc)
1 static char help[] = "Test PETSc integer hash map.\n\n";
2 
3 #include <petsc/private/hashmapi.h>
4 #include <petsc/private/hashmapiv.h>
5 #include <petscsys.h>
6 
7 /* Unused, keep it for testing purposes */
8 PETSC_HASH_MAP(HMapIP, PetscInt, void*, PetscHashInt, PetscHashEqual, NULL)
9 
10 /* Unused, keep it for testing purposes */
11 typedef struct { double x; double y; double z; } Point;
12 static Point origin = {0.0, 0.0, 0.0};
13 PETSC_HASH_MAP(HMapIS, PetscInt, Point, PetscHashInt, PetscHashEqual, origin)
14 
15 #define PetscAssert(expr) do {            \
16 if (PetscUnlikely(!(expr)))               \
17   SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB, \
18            "Assertion: `%s' failed.",     \
19            PetscStringize(expr));         \
20 } while (0)
21 
22 int main(int argc,char **argv)
23 {
24   PetscHMapI     ht = NULL, hd;
25   PetscHMapIV    htv;
26   PetscInt       n, v, koff, keys[4], voff, vals[4],na,nb,i,size,*karray,off;
27   PetscScalar    *varray,*vwork;
28   PetscBool      has, flag;
29   PetscErrorCode ierr;
30 
31   ierr = PetscInitialize(&argc,&argv,NULL,help);if (ierr) return ierr;
32 
33   ierr = PetscHMapICreate(&ht);CHKERRQ(ierr);
34   PetscAssert(ht != NULL);
35   ierr = PetscHMapIGetSize(ht,&n);CHKERRQ(ierr);
36   PetscAssert(n == 0);
37 
38   ierr = PetscHMapIResize(ht,0);CHKERRQ(ierr);
39   ierr = PetscHMapIGetSize(ht,&n);CHKERRQ(ierr);
40   PetscAssert(n == 0);
41 
42   ierr = PetscHMapIHas(ht,123,&has);CHKERRQ(ierr);
43   PetscAssert(has == PETSC_FALSE);
44   ierr = PetscHMapIGet(ht,123,&v);CHKERRQ(ierr);
45   PetscAssert(v == -1);
46 
47   ierr = PetscHMapISet(ht,123,42);CHKERRQ(ierr);
48   ierr = PetscHMapIGetSize(ht,&n);CHKERRQ(ierr);
49   PetscAssert(n == 1);
50   ierr = PetscHMapIHas(ht,123,&has);CHKERRQ(ierr);
51   PetscAssert(has == PETSC_TRUE);
52   ierr = PetscHMapIGet(ht,123,&v);CHKERRQ(ierr);
53   PetscAssert(v == 42);
54 
55   ierr = PetscHMapIDel(ht,123);CHKERRQ(ierr);
56   ierr = PetscHMapIGetSize(ht,&n);CHKERRQ(ierr);
57   PetscAssert(n == 0);
58   ierr = PetscHMapIHas(ht,123,&has);CHKERRQ(ierr);
59   PetscAssert(has == PETSC_FALSE);
60   ierr = PetscHMapIGet(ht,123,&v);CHKERRQ(ierr);
61   PetscAssert(v == -1);
62 
63   ierr = PetscHMapIQuerySet(ht,123,1,&flag);CHKERRQ(ierr);
64   PetscAssert(flag == PETSC_TRUE);
65   ierr = PetscHMapIQuerySet(ht,123,1,&flag);CHKERRQ(ierr);
66   PetscAssert(flag == PETSC_FALSE);
67   ierr = PetscHMapIQueryDel(ht,123,&flag);CHKERRQ(ierr);
68   PetscAssert(flag == PETSC_TRUE);
69   ierr = PetscHMapIQueryDel(ht,123,&flag);CHKERRQ(ierr);
70   PetscAssert(flag == PETSC_FALSE);
71 
72   ierr = PetscHMapIResize(ht,13);CHKERRQ(ierr);
73   ierr = PetscHMapIGetSize(ht,&n);CHKERRQ(ierr);
74   PetscAssert(n == 0);
75 
76   ierr = PetscHMapIClear(ht);CHKERRQ(ierr);
77   ierr = PetscHMapIGetSize(ht,&n);CHKERRQ(ierr);
78   PetscAssert(n == 0);
79 
80   ierr = PetscHMapISet(ht,321,24);CHKERRQ(ierr);
81   ierr = PetscHMapISet(ht,123,42);CHKERRQ(ierr);
82   ierr = PetscHMapIGetSize(ht,&n);CHKERRQ(ierr);
83   PetscAssert(n == 2);
84 
85   koff = 0; keys[0] = keys[1] = 0;CHKERRQ(ierr);
86   ierr = PetscHMapIGetKeys(ht,&koff,keys);CHKERRQ(ierr);
87   ierr = PetscSortInt(koff,keys);CHKERRQ(ierr);
88   PetscAssert(koff == 2);
89   PetscAssert(keys[0] == 123);
90   PetscAssert(keys[1] == 321);
91 
92   voff = 0; vals[0] = vals[1] = 0;
93   ierr = PetscHMapIGetVals(ht,&voff,vals);CHKERRQ(ierr);
94   ierr = PetscSortInt(voff,vals);CHKERRQ(ierr);
95   PetscAssert(voff == 2);
96   PetscAssert(vals[0] == 24);
97   PetscAssert(vals[1] == 42);
98 
99   koff = 0; keys[0] = keys[1] = 0;
100   voff = 0; vals[0] = vals[1] = 0;
101   ierr = PetscHMapIDuplicate(ht,&hd);CHKERRQ(ierr);
102   ierr = PetscHMapIGetKeys(ht,&koff,keys);CHKERRQ(ierr);
103   ierr = PetscHMapIGetVals(ht,&voff,vals);CHKERRQ(ierr);
104   ierr = PetscSortInt(koff,keys);CHKERRQ(ierr);
105   ierr = PetscSortInt(voff,vals);CHKERRQ(ierr);
106   PetscAssert(koff == 2);
107   PetscAssert(voff == 2);
108   PetscAssert(keys[0] == 123);
109   PetscAssert(keys[1] == 321);
110   PetscAssert(vals[0] == 24);
111   PetscAssert(vals[1] == 42);
112   ierr = PetscHMapIDestroy(&hd);CHKERRQ(ierr);
113 
114   ierr = PetscHMapISet(ht,0,0);CHKERRQ(ierr);
115   ierr = PetscHMapIGetSize(ht,&n);CHKERRQ(ierr);
116   PetscAssert(n != 0);
117   ierr = PetscHMapIReset(ht);CHKERRQ(ierr);
118   ierr = PetscHMapIGetSize(ht,&n);CHKERRQ(ierr);
119   PetscAssert(n == 0);
120   ierr = PetscHMapIReset(ht);CHKERRQ(ierr);
121   ierr = PetscHMapIGetSize(ht,&n);CHKERRQ(ierr);
122   PetscAssert(n == 0);
123   ierr = PetscHMapISet(ht,0,0);CHKERRQ(ierr);
124   ierr = PetscHMapIGetSize(ht,&n);CHKERRQ(ierr);
125   PetscAssert(n != 0);
126 
127   ierr = PetscHMapIDestroy(&ht);CHKERRQ(ierr);
128   PetscAssert(ht == NULL);
129 
130   ierr = PetscHMapICreate(&ht);CHKERRQ(ierr);
131   ierr = PetscHMapIReset(ht);CHKERRQ(ierr);
132   ierr = PetscHMapIGetSize(ht,&n);CHKERRQ(ierr);
133   PetscAssert(n == 0);
134   ierr = PetscHMapIDestroy(&ht);CHKERRQ(ierr);
135 
136   ierr = PetscHMapIVCreate(&htv);CHKERRQ(ierr);
137   n = 10;
138   ierr = PetscHMapIVResize(htv,n);CHKERRQ(ierr);
139   ierr = PetscHMapIVGetCapacity(htv,&na);CHKERRQ(ierr);
140   PetscAssert(na>=n);
141   for (i=0; i<n; i++) {
142     ierr = PetscHMapIVSet(htv,i+100,10.);CHKERRQ(ierr);
143   }
144   ierr = PetscHMapIVGetCapacity(htv,&nb);CHKERRQ(ierr);
145   PetscAssert(nb>=na);
146   for (i=0; i<(2*n); i++) {
147     ierr = PetscHMapIVAddValue(htv,i+100,5.);CHKERRQ(ierr);
148   }
149   ierr = PetscHMapIVGetSize(htv,&size);CHKERRQ(ierr);
150   PetscAssert(size==(2*n));
151   ierr = PetscMalloc3(size,&karray,size,&varray,size,&vwork);CHKERRQ(ierr);
152   off = 0;
153   ierr = PetscHMapIVGetPairs(htv,&off,karray,varray);CHKERRQ(ierr);
154   PetscAssert(off==(2*n));
155   ierr = PetscSortIntWithDataArray(off,karray,varray,sizeof(PetscScalar),vwork);CHKERRQ(ierr);
156   for (i=0; i<n; i++) {
157     PetscAssert(karray[i]==(i+100));
158     PetscAssert(karray[n+i]==(n+i+100));
159     PetscAssert(varray[i]==15.);
160     PetscAssert(varray[n+i]==5.);
161   }
162   ierr = PetscFree3(karray,varray,vwork);CHKERRQ(ierr);
163   ierr = PetscHMapIVDestroy(&htv);CHKERRQ(ierr);
164 
165   ierr = PetscFinalize();
166   return ierr;
167 }
168 
169 /*TEST
170 
171    test:
172 
173 TEST*/
174